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Abstract — We consider the probability for a colour-singlet qq pair to emit a gluon, in strongly 
and smoothly ordered antenna showers. We expand to second order in Ug and compare to the 
second-order QCD matrix elements for Z — > 3 jets, neglecting terms suppressed by 
We give a prescription that corrects the shower to the matrix-element result at this order for 
both soft and hard emissions, thereby explicitly reducing its dependence on evolution- and 
renormalization-scale choices. We confirm that the choice of p± for both of these scales ab- 
sorbs all logarithms through 0{a'^), and contrast this with various alternatives. We include 
these corrections in the VINCIA shower generator and study the impact on LEP event-shape 
and fragmentation observables. An uncertainty estimate is provided for each event, in the form 
of a vector of alternative weights. 



1 Introduction 

Modern QCD descriptions of hard-scattering events at particle colliders can be roughly divided into two 
broad categories. In the first, fixed-order descriptions, matrix elements are computed for all allowed 
initial states with a given final state, F, plus a limited number of additional partons. The leading-order 
(LO) description has the minimal number (often zero) of additional partons. For improved accuracy, 
one includes matrix elements with one extra parton beyond leading order and one loop correction (next- 
to-leading order) and so forth. The squared matrix elements are numerically integrated over the allowed 
phase space, after accounting for any divergences. Given the accuracy, i.e. order of the description, 
the possible number of additional final-state particles is in essence predetermined, and can take one 
(LO), or two (NLO) etc, values. In the second, parton-shower descriptions, one also starts from matrix 
elements for the desired hard process, F, but additional radiation is now generated stochastically via a 
shower algorithm, which is essentially Markovian. This is a unitary process, with probability one, and 
therefore does not change the probability of the underlying hard process to occur. The number of final- 
state partons is now not predetermined, and can take an infinity of different values. The two approaches 
have complementary strengths and weaknesses (for a review, see e.g. [1]). When hard extra emissions 
(e.g., hard jets) are important to model well, one looks to descriptions of the first category. However, 
the calculation is then unpredictive for near-collinear and soft radiation (e.g., jet substructure and soft 
wide-angle jets). The obverse holds for the second category. 

Even from this very cursory summary, it is clear that methods to unite the two — combining 
strengths and eliminating weaknesses — are very important. Two longstanding and very successful 
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approaches for combining one-loop matrix elements with parton showers are MC@NLO [2-6] and 
POWHEG [7-9]. An important restriction of both of these is that only the spectra of the basic hard- 
scattering partons are corrected to NLO precision, while those of additional QCD emissions are not. 
Removing this limitation, fully or partially, has been the focus of much recent effort [10-13], and is also 
among the main goals of this paper. While most other approaches employ parton-shower algorithms 
which are based on 1 — )• 2 splitting kernels, we develop an approach for matching NLO descriptions 
to showers based on 2 — 3 splittings [14]. The equivalent of the 1 — 2 splitting kernels are, for us, 
dipole-antenna functions [14-16]. At the practical level, our approach is in the context of the VINCIA 
framework [17, 18]. Whatever the splitting kernel, the parton-shower approaches rest upon the factor- 
ization of both phase space and matrix element when the splitting is either soft or coUinear, or both. 
A technical advantage of our approach is that the (n + 3) -particle phase-space factorizes exactly into a 
(n + 2) -particle phase space times a 2 — 3 phase space with all momenta on-shell, without need for mo- 
mentum reshuffling [19]. Phase-space factorization and the antenna-based matrix-element factorization 
are important to our approach in about equal measure. 

The essential bottleneck in such combinations of fixed order and parton shower is how to avoid 
double counting both real emissions as well as virtual effects. A key aspect of this is how well the NLO 
emissions are mimicked by parton-shower emissions. Emissions generated by a parton-shower Markov 
chain in fact produce approximations to tree-level matrix elements up to arbitrary numbers of legs, while 
the no-emission Sudakov factors generate the equivalent all-orders loop corrections'. This all-orders 
resummation of contributions is ordered in a measure of jet resolution, called the evolution scale, which 
we denote Qe- It is typically chosen to be a measure of transverse momentum [14] or invariant mass. 
Its fundamental role is to separate resolved from unresolved emissions, in analogy to a jet-clustering 
measure. The different evolution variables each have their strengths, depending on the context. As part 
of our present study, we judge them by how well their fixed-order expansions approximate the NLO 
matrix elements. 

The main purpose of this paper is to define, for e+e" initial states, an antenna-based shower al- 
gorithm that incorporates multileg NLO corrections for both soft and hard emissions, and to study the 
quality of the matching for a variety of evolution variables. We strive for next-to-leading logarithmic 
(NLL) accuracy, in a way we shall detail below. 

The leading-logarithmic (LL) structure of antenna showers was discussed in [21,22], with explicit 
comparisons of various algorithms to tree-level 0(0^) matrix elements presented in [18,22,23]. A pre- 
scription for matching the showers to reproduce tree-level matrix elements exactly (over all of phase 
space) was developed in [18], with uncertainty variations provided in the form of vectors of alternative 
weights for each event. In [24] and [25] substantial speedups of the matching algorithm were obtained 
by dividing phase space into so-called sectors, and by deriving a formalism for using individual helic- 
ity amplitudes to correct the shower evolution, respectively. To further probe the structure of antenna 
showers, at the subleading-logarithmic level, we shall here consider the expansion of exclusive 2 — )• 3 
splitting probabilities to ©(a^), comparing these to one-loop matrix elements [26] and to corresponding 
second-order antenna functions [16,27]. 

We shall compare six different types of ordering criteria for the shower evolution: 1) strong ordering 
in transverse momentum, 2) strong ordering in dipole virtuality, 3 & 4) strong ordering in two variants 
of emission energy (mostly intended as cross-checks), and 5 & 6) so-called smooth ordering in p± and 

'For an introduction to such ciiains and a description of tlieir properties, see e.g. [1, 20]. 
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in dipole virtuality, as defined by [18]. We also consider several different choices for the renormalization 
scale ij.f( used in the tree-level antenna functions and discuss how to systematically absorb contributions 
proportional to the /3-function by this choice, elaborating on earlier arguments [28, 29]. 

Finally, we will present a prescription for how to systematically incorporate the second-order (one- 
loop) qq — )• qgq antenna into the shower evolution, for each of the studied evolution variable choices. 
This will essentially constitute the NLL accuracy mentioned above. The resulting shower algorithm, 
whose qq — qgq splitting probability should therefore be correct to 0{ag) over all of phase space, has 
been implemented in the publicly available VINCIA plug-in [17] to the PYTHIA 8 event generator [30]. 

We have organized the paper as follows. In section 2 we discuss introductory aspects of (antenna) 
shower algorithms, define the various evolution variables, and the implementation of an ordering pre- 
scription that rules the shower evolution. In section 3 we present our matching prescription in detail, 
initially for 3-parton final states in Z-decay, then generally for n partons. In section 4 we discuss details 
of the Sudakov integrals required in the matching prescription and compare the infrared hmits of those 
integrals to those of the one-loop matrix elements. In section 5 we combine one-loop and tree-level cor- 
rections in a single algorithm, perform speed benchmarking, and study the impact on LEP observables, 
especially in the context of as{mz) extractions. We conclude in section 6 and elaborate on technical 
aspects in the appendices. 

2 Antenna Showers 

In this section, we recap the basic antenna-shower formalism, as used in the viNCiA shower algorithm. 
This also serves to introduce the basic notation and conventions that will be used in later sections. 

2.1 The Formal Basis of Antenna Showers 

Antenna showers are based on the factorization of (squared) colour-ordered QCD amplitudes in soft and 
collinear limits, which can be expressed as follows 

|M(. . . ,p,,p„ . . 'il^ C ^ |M(. . . ,K +p„ . . (1) 

Sij 

\M{...,pi,pj,pk,...)\ g^C Ag{sij,Sjk,Sijk)\M{...,pi,pk,...)\ , (2) 

with g"^ = AttUs the strong coupling and the subscript g in the second line emphasizing that the soft 
Umit is only relevant for gluons. 

In the collinear limit (first line), P{z) are the Altarelli-Parisi splitting kernels [31], z is the energy 
fraction taken by parton i (with a fraction (1 — z) going to parton j), and C is a colour factor, which we 
discuss below. This hmit forms the basis for traditional parton showers, such as those in the PYTHIA 
generator [32]. 

In the soft-gluon limit (second line), the function A has dimension GeV~^, and is called an antenna 
function. For unpolarized massless partons^, its leading term is the so-called eikonal or dipole factor, 

AEik{Sij , Sjk, Sijk) = '- , (3) 
"In the context of massive particles, replace Sat by 2pa ■ pt in all expressions. For a more complete treatment, see [33]. 
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where Sik = Sijk — Sij — Sjk for massless partons. It was found early on that this factor can be repro- 
duced by a traditional parton shower by imposing the requirement of angular ordering of subsequent 
emissions [34]. This gave rise to the angular-ordered showers [35,36] in the HERWIG and HERWlG-i-i- 
generators [37,38] as well as the imposition of an angular-ordering constraint [32,39] in the JETSET and 
PYTHIA generators [30,40]. 

In fixed-order calculations, dipole [41] and antenna [15, 16, 19] functions are frequently used to 
define subtraction terms. These functions include additional subleading terms, beyond the eikonal one, 
which are necessary to correctly describe both soft and coUinear limits in all regions of phase space. In 
the parametrization we shall use, their most general forms, for the branching process IK — ijk, are 

Aspiit(,Sij,Sjk,m'jj^) = — 2~ i ^ + -fspiit) , (5) 

for gluon-emission and gluon-spUtting processes, respectively, with the parent antenna invariant mass, 
"^/E" — {Pi ^Pk)"^ = {Pi +Pj + PA:)^ and the scaled invariants, 

^ij _ ^jk 

Vij — 9 ) Uik — 9 ) \P) 
^IK ^IK 

and we use the notation 6ig = 1 if parton i is a gluon and zero otherwise. The functions FEmit and Fspiit 
allow for the presence of non-singular terms, which are in principle arbitrary. A logical choice would 
be F = 0, but this would not be invariant under reparametrizations of the antenna functions across the 
gluon-coUinear singular limits [18]. Since the F functions can anyway be made useful in the context of 
uncertainty estimates [17, 18], we therefore leave them as functions whose forms we are free to choose. 

In the soft-gluon limit, the eikonal factor is reproduced by the first term in in eq. (4). In the collinear 
q ^ qg limit, the AP splitting kernel also is reproduced. For collinear g ^ gg and g ^ qq branchings, 
one must sum over the contributions from two neigbouring antennae, which together reproduce the AP 
splitting kernel. Limits that are both soft and collinear are also correctly reproduced [16]. 

In the antenna context, the colour factors are 2Cf for qq — )• qgq, Ca for gg — )• ggg^, and 2Tji for 
gluon splitting to qq, again using the normalization convention adopted in [18]. However, for qg — qgg 
there is a genuine subleading ambiguity whether to prefer, say, 2Cf, Ca, or something interpolating 
between them [42] . At fixed order, the question of subleading colour can in fact be dealt with quite ele- 
gantly, by using Ca for all antennae and then including an additional qq antenna with a negative colour 
factor, —Ca/Nq, spanned between the two endpoint quarks, for each qg . . .q chain [43]. In the context 
of an antenna-based shower, however, it is desirable to use only positive-definite antenna functions, and 
a prescription for absorbing the negative one into the positive ones was given in [18]. In the context of 
this work, however, we shall largely ignore subleading-colour aspects and, unless explicitly stated oth- 
erwise, assign a colour factor Ca to the qg — )• qgg antenna function, thereby overcounting the collinear 
limit in the quark direction by a factor Ca/ {2Cf) ^ 1 + ^/N^. 



^^Note that in a process like — >■ gg, there are two gg antennae at the Born level, and hence the antenna approximation to 
H° — >■ ggg is twice as large as the single gg — >■ ggg antenna. Likewise, the collinear limit of a gluon is obtained by summing 
over the contributions from both of the dipoles/antennae it is connected to. One must also include a sum over permutations of 
the final-state gluons, if comparing to a summed matrix-element expression. 
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The renormalization scale used to evaluate the strong coupling in the antenna function, g1 = 
47ras(/ips), is typically chosen proportional to (following [28]). As alternatives, we shall also con- 
sider /Xpg oc mj-) = 2min(sij, Sjk), and, as an extreme case which connects with fixed-order calcula- 
tions, the invariant mass of the antenna, ;Upg tx m?jj^. 

A final aspect concerns the phase-space factorization away from the collinear limit. Within the 
framework of collinear factorization (and hence, in traditional parton showers), the momentum fraction, 
z, is only uniquely defined in the exactly collinear limit; outside that limit, the choice of z is not unique. 
In addition, a prescription must be adopted for ensuring overall momentum conservation, leading to the 
well-known ambiguities concerning recoil strategies (see e.g. [1]). In antenna showers, on the other 
hand, the antenna function is defined in terms of the unique branching invariants, sij and Sjk, over all 
of phase space, and the phase space itself has an exact Lorentz-invariant and momentum-conserving 
factorization, 

= d$n-l X d^>Ant , (7) 

with 

1 

= d..,d.,,- (8) 

for massless partons (for massive ones, see [33]), with the (j) angle parametrizing rotations around the 
antenna axis, in the CM of the antenna. Note the equality signs; no approximation is involved in this 
step. The only remaining phase-space ambiguity, outside the singular limits, is present when specifying 
how the post-branching momenta are related to the pre -branching ones. This is defined by a kinematics 
map, the antenna equivalent of a recoil strategy, which we here take to be of the class defined by [15, 17]. 

2.2 Constructing a Shower Algorithm 

In a shower context, the amplitude and phase-space factorizations above imply that we can interpret the 
radiation functions (AP splitting kernels or dipole/antenna functions) as the probability for a radiator 
(parton or dipole/antenna) to undergo a branching, per unit phase-space volume, 

^=.?C^W, (9) 

where we use <I> as shorthand to denote a phase-space point. (If there are several partons/dipoles/antennae, 
the total probability for branching of the event as a whole is obtained as a sum of such terms.) 

An equally fundamental object in both analytical resummations and in parton showers is the Sudakov 
form factor, which defines the probability for a radiator not to emit anything, as a function of the shower 
evolution parameter, Q (i.e., similarly to a jet veto, with Q playing the role of the jet clustering scale; 
we return to the choice of functional form for the shower evolution scale in section 2.3). In the all- 
orders shower construction, these factors generate the sum over virtual amplitudes plus unresolved real 
radiation, and hence their first-order expansions play a crucial role in matching to next-to-leading order 
matrix elements. We here recap some basic properties. The Sudakov factor, giving the no-emission 
probability between two values of the shower evolution parameter, Qi and Q2 (with Qi > Q2), is 
defined by 

A{Ql Ql) = exp (- j^^ d$ j = exp (- j^^ gl C A{^) d$ j , (10) 
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where it is understood that the integral boundaries must be imposed either as step functions on the 
integrand or by a suitable transformation of integration variables, accompanied by Jacobian factors. 

This description has a very close analogue in the simple process of nuclear decay, in which the 
probability for a nucleus to undergo a decay, per unit time, is given by the nuclear decay constant, 

-^ = CN. (11) 
The probability for a nucleus existing at time ti to remain undecayed before time t2, is 



A(t2, ii) = exp ^— y cn dt^ = exp(— CivAt) . 



(12) 



This case is especially simple, since the decay probability per unit time, cat, is constant. By conservation 
of the total number of nuclei (unitarity), the activity per nucleon at time t, equivalent to the "resummed" 
decay probability per unit time, is minus the derivative of A, 

dPresjt) dA 

In QCD, the emission probability varies over phase space, hence the probability for an antenna not to 
emit has the more elaborate integral form of eq. (10). By unitarity, the resummed branching probability 
is again minus the derivative of the Sudakov factor, 

= g',CA{<^>)A{QlQ\'^)), (14) 

with the shower evolution scale (typically chosen as a measure of invariant mass or transverse 

momentum, see section 2.3), evaluated at the phase-space point 

In shower algorithms, branchings are generated with this distribution, starting from a uniformly 
distributed random number TZ G [0, 1], by solving the equation, 

TZ = A{QlQ^) (15) 

for Q^. For an initial distribution of "trial" branching scales, we do not use the full antenna function, 
eq. (4), as the evolution kernel, but only its leading singularity, 

A. = 5!!k = ^, (16) 

SijSjk 

where p±A is the ARIADNE definition of transverse momentum [44], which is also the one used in 
VINCIA. This reflects the universal l/p^ behaviour of soft-gluon emissions. In addition to the trial 
scale, Q, two complementary phase-space variables are also generated (which we usually label ( and 
[18]), according to the shape of At over a phase-space contour of constant Q. From these, the model- 
independent set of trial phase-space variables (sij, sjk, (p) are determined by inversion of the defining 
relations Q{sij, Sjk) and C(Sij) Sjfc)> and the full kinematics (i.e., four-momenta) of the trial branching 
can then be constructed [17]. 
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To decide whether to accept the trial or not, we note that the function At differs from the eikonal in 
eq. (3) by the replacement of Sik in the numerator by m'jj^. By accepting the trial scales generated by 
At with the probability 

Peik = ^ = ^ < 1 , (17) 

At mjj^ 

the eikonal approximation can be recovered, by virtue of the veto algorithm [1,45,46]. Crucially, any 
other function that has the eikonal as its soft-collinear limit could equally well be imprinted on the trial 
distribution by a similar veto. Two particularly relevant choices are the full physical antenna function, 
eq. (4) (which includes additional coUinear-singular terms in addition to the eikonal) and the GKS- 
corrected antenna function (which also incorporates a multiplicative factor that represents tree-level 
matching in VINCIA), 

pLL _ ^Emit ,^oN 

-'accept ) \^"J 

pLO _ pLL T-> /"I Q\ 

-'accept -'accept -'^n ; \^-'J 

with ^Emit given in eq. (4) and Rn the n-parton tree-level GKS matching factor [18], to which we return 
in section 3.1. 

Note that, for gluon-splitting antenna functions (Xg — Xqq), we use Q = niqq, with a trial function 
oc l/m^g, and again implement the physical antenna function, eq. (5), and LO matching corrections by 
vetos. We also include the so-called ariadne factor, PAri> defined by 

2S7V 

^Split — ^ -PAri ^Split = ; ^Split , (20) 

sp + Sn 

with Sn the invariant mass squared of the colour neighbour on the other side of the splitting gluon 
and Sp = m'jj^ the invariant mass squared of the parent (splitting) antenna. This does not modify 
the singular behavior (as will be elaborated upon below), and was shown to give significantly better 
agreement with the Z — qqqq matrix element in [24]. 

Explicit solutions to eq. (15) using the trial function defined by eq. (16) were presented in [18], for 
fixed and first-order running couplings. In the context of the present work, two-loop running has been 
implemented using a simple numerical trick, as follows: given a value of as{Mz), we determine the 
corresponding two-loop value of Aq"^^^"^. We then use that A value in the one-loop solutions in [18], 
and correct the resulting distribution by inserting an additional trial accept veto: 

a^'°°P(Q,A?^'°°P) 
p2— loop cjcu / 

-'"^accept — 1-loop/^ . 2-loopN ' 

as p(Q,AqcdP) 

Due to the faster pace of 2-loop running, a'i^^°°^{Q, A) < al^^°°^{Q, A), hence this accept probability 
is guaranteed to be smaller than or equal to unity. 

A final point concerns if there are several "competing" radiators (equivalent to several competing 
nuclei, and/or several competing available decay channels for each nucleus). In this case, the trial 
with the highest value of Q is selected (corresponding to the one happening at the earliest time, t), 
and consideration of any other branchings (decays) are postponed temporarily. After a branching, any 
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partons involved in that branching are replaced by the post-branching ones, and any postponed trial 
branchings involving those partons are deleted. The evolution is then restarted, from the scale Q of the 
new configuration, until there are no radiators left with trial branching scales larger than a fixed, lower, 
cutoff, normally identified with the hadronization scale, Qhad ~ 1 GeV. 



2.3 Evolution and Ordering 

In order to solve eq. (15) we need to specify the form of eq. (10), which takes us from one scale Q\ 
to a lower scale Q^. We change variables to parametrize the integral by the ordering variable, Q, and 
another, complementary (but otherwise arbitrary), phase-space variable which we denote by C,. The 
generic evolution integral now reads 

rQ\ 

A{QlQl) = / ^ glCdQ^dC\J\A{Q\0 (22) 

with I J| denoting the Jacobian of this transformation. For branchings involving gluon emission, we con- 
sider three possible choices for the ordering variable: dipole virtuality niD, transverse momentum, and 
the energy of the emitted parton, E* (in the CM of the parent antenna), with the following definitions. 



Qei 


= ml) 






(23) 


Q'e2 


= 4pl 


= 4mf^ yijyjk , 




(24) 


Qe3 


= AEf 


= iVij + Vjkf 


2 2 

— Xj rrijj^ , 


(25) 



with the energy fraction Xj = 2E* /mix- 

All three options are available as ordering variables in the VINCIA shower Monte Carlo. They are 
illustrated in fig. 1, where contours of constant value ofyE = Ql/fn^^ are shown for each variable, as 
a function of yij and yj^. For completeness, we show both the case of a linear (top row) and quadratic 
(bottom row) dependence on the branching invariants, for each variable. Since the ordering variable 
raised to any positive power will result in the same relative ordering of emissions within a given antenna, 
the distinction between linear and quadratic forms does not affect individual antenna Sudakov factors. 
It does, however, affect the "competition" between different antennae, and the choice of restart scale for 
subsequent evolution after a branching has taken place, as will be discussed further below. 

In labeling the columns in fig. 1, we have also emphasized that mass-ordering, as defined here, corre- 
sponds to choosing the smallest of the daughter antenna masses as the "resolution scale" of the branch- 
ing, whereas p_\_ and energy correspond to using the geometric and arithmetic means of the daughter 
invariants, respectively. Naively, each of these could be taken as a plausible measure of the resolution 
scale of a given phase-space point. We shall see below which ones lead to better agreement with the 
one-loop matrix elements. 

We consider two possible definitions for the complementary phase-space variable C,, 

Ci = ^2 (26) 

Vij + Vjk 

C,2 = Vij ■ (27) 
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Figure 1 : Contours of constant value of y^; = Q'^/m'jj^ for evolution variables linear ( top) and quadratic 
(bottom) in the branching invariants, for virtuality-ordering (left), -ordering (middle), and energy- 
ordering (right). Note that the energy-ordering variables intersect the phase-space boundaries, where 
the antenna functions are singular, for finite values of the evolution variable. They can therefore only be 
used as evolution variables together with a separate infrared regulator, such as a cut in invariant mass, 
not shown here. 



We emphasize that the choice of ( has no physical consequences, it merely serves to reparametrize the 
Lorentz-invariant phase space. We may therefore let the choice be governed purely by convenience, and, 
for each antenna integral, select whichever of the above definitions give the simplest final expressions. 
The corresponding Jacobian factors, for each of the evolution-variable choices we shall consider, are 
listed in tab. 1. 

Note that, for the special case of the mj-, and mjj variables, which contain the non-analytic function 
min(yjj, yjk), the C definitions in eqs. (26) and (27) apply to the branch with yij > yjk- For the other 
branch, yij and yjj^ should be interchanged. With this substitution, the Jacobians listed in tab. 1 apply to 
both branches'*. 

For branchings involving gluon splitting, g — )• qq, we restrict our attention to two possibilities, 
''This corresponds to replacing yij by max(yij, yjk) in the numerator of eq. (26) and in eq. (27). 
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Note: \JiQ\C)\ = Sijk\J{yE,0\ 



Table 1: Jacobian factors for all combinations of evolution variables and ( choices. 

ordering in p±, defined as above, and ordering in gluon virtuality, defined as 

Qe4 = m^g* = nigg (for gluon splitting) . (28) 

Note that the normalization factors for the ordering variables have in all cases been chosen such that the 
maximum value of the ordering variable is mjj^. 

Since the phase space for subsequent branchings is limited both by the scale Qe of the previous 
branching (according to strong ordering) and by the invariant mass of the antenna rrij, the effective 
"restart scale", after a branching in a strongly ordered shower, is given by 

Q% = mm{Q\mJ) , (29) 

for each antenna j. 

Depending on the choice and value of Q, one or both daughter antennae after a splitting may have 
a non-trivial restriction on the phase space available for subsequent branching. Conversely, if Q > mj, 
there is no such restriction. Physically, we distinguish between the case in which the strong-ordering 
condition implies a non-trivial constraint on the evolution of the produced antennae, eating into the 
phase-space that would otherwise be accessible, and the case in which the strong-ordering condition 
does not imply such a constraint. 

The regions of qq — qgq phase space in which either zero, one, or both of the daughter antennae 
{qg and gq respectively) are constrained by the ordering condition are illustrated in fig. 2, for each of 
the choices of evolution variable under consideration. The black shaded areas correspond to regions in 
which both the qg and gq antennae are restricted, by having Q < nij. The darker gray shaded areas show 
regions in which only one of the antennae is restricted, while the other will still be allowed to evolve 
over its full phase space. In the light-gray shaded areas, both of the antennae are allowed to evolve 
over all of their available phase spaces, equivalent to the ordering condition imposing no constraint on 
the subsequent evolution. We recall that we are here discussing the upper boundary on the subsequent 
evolution, hence the infrared^ (IR) poles are not affected. 

To further clarify the meaning of the plots in fig. 2, let us discuss panel (e) as an example. The 
coordinates, {yij,yjk), represent the 3-parton state before it evolves to a 4-parton state, and each point 
corresponds to a specific value of the evolution variable at hand, cf. fig. 1. Assuming ordering in p± and 
using subscript (3) for quantities evaluated in the 3-parton state, the value of the evolution variable for a 
specific {uij , yj^) point is (5|;(3) = 4p5_(3) = ^yijyjkS, with s = m\ at the Z pole. The further evolution 
of the shower, from a 3- to a 4-parton state, involves a sum over all possible branchings of the qg and 

'Note: we use the word infrared to refer collectively to soft and/or coUinear regions of phase space. 
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Figure 2: Illustration of the regions of 3-parton phase space in which the subsequent evolution of the 
qg and gq antennae is restricted (from above) by the strong-ordering condition. See the text for further 
clarification of this plot. Black: both antennae restricted. Dark Gray: one antenna restricted, the other 
unrestricted. Light Gray: both antennae unrestricted. Top/Bottom: linear/quadratic in the branching 
invariants, for mass-ordering (left), -ordering (middle), and energy-ordering (right). 



gq antennae. Consider the qg one. Its branchings can again be characterized by two invariants (si, S2), 
both of which will be smaller than rrigg. However, depending on the value of rn^g (or, equivalently, i/ij) 
the p± of the new configuration, 4p^^^^ = AsiS2/m'^g may actually be larger than 4:p']_(^^y In a strongly 
ordered shower, such configurations are not allowed, and would be discarded. Whether this situation 
can occur or not, for one or both of the qg and gq antennae, as a function of {yij,yjk), is what fig. 2 
reveals, for each type of ordering variable. 

The mathematical consequence is that, in the dark- and black-shaded regions, respectively, the upper 
boundary of one or both of the qg and gq antenna integrals is set by the evolution variable, rather than by 
phase space. This creates an important difference between the integrals generated by a shower algorithm 
and those used for IR subtractions in traditional fixed-order applications (for which the integrals always 
run over all of phase space). In particular, we see that the strong-ordering condition will generate 
additional logarithms involving •Sij /<5|;(3) argument. For a "good" choice of evolution variable, these 
logarithms should explicitly cancel against ones present in the one-loop matrix elements, a question we 
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Ordering type 


Ql 


CramiQl) Cmax(Q|;) 


3^4 restriction 


, . linear 
p_L -ordering 

squared 

, . linear 
-ordering 

squared 

rn* J • linear 
i? -ordering 

squared 


^m]j^ yijVjk 

4m2^ min(?/2.,y2^) 

iVij + yjk) 


1 V ^ h, ' IK 

2 

2 ^ 2 

/ Qe ^ / Qe 

V "'"^K V ""''^ 
1 

1 


& {ml^t -2min{sij,Sjk)) 

^ (^mL - 4 j 
1 



Table 2: Boundaries corresponding to the ordering variables portrayed in fig. 1, with ml^^^ corresponding 
to the active 3 — )• 4 dipole Sqg or Sgq, and s = m| at the Z pole. We have chosen ^2 as the energy 
sharing variable for rriD and p± ordering and (i for E* ordering, with ( defined as in eq. (26) and 
eq. (27). The energy variable will lead to infinities if the hadronization scale is not imposed as a cut-off. 



shall return to in detail in section 4. 

Several interesting structures can be seen in fig. 2. Firstly, the linearized variables imply less severe 
constraints on the subsequent evolution than the quadratic ones. This is easy to understand given that 
the linearized variables, Qnn, are related to the quadratic ones, Qqdr. by 

Qlin = Qqdv miK , (30) 

and hence Qun > <5qdr> implying a higher absolute restart scale for the linearized ordering variables. 

It is also apparent that, for a given choice of linearity, mass-ordering reduces the phase-space for 
further evolution more than -ordering does, which in turn is more constraining than energy-ordering. 
In this comparison, however, it becomes important to recall that the traditional ordering variables used, 
e.g., in VINCIA, are the linearized mass-ordering and the quadratic p±_ and energy-ordering variables^. 
Within that group, p^-ordering appears to be the most restrictive, followed by energy-ordering, with 
(traditional, linearized) mass-ordering leading to the most open phase space for the subsequent evolu- 
tion. 

We are now able to fully specify the boundaries of the evolution integrals in eq. (22). For each Qe 
contour (see fig. 1), the integration limits in ( are listed in tab. 2. Combined with slQe interval and an 
antenna function, these boundaries account for the integrated tree-level splitting probability when going 
from one scale to another Q2, as expressed by eq. (22). The last column in tab. 2 tells when the 
3 — >^ 4 ordering condition is active. In fig. 2 this corresponds to a region darkening due to the restriction, 
with its shade determined by the amount of restricted dipoles. 

Finally, we note that the dependence on Q in eq. (29) causes explicit non-Markovian behavior at the 
4-parton level and beyond, since the value of Q then depends explicitly on which branching was the last 
to occur. A more strictly Markovian variant of this is obtained by letting the min() function act on all 
possible Q values (corresponding to all possible colour-connected clusterings) of the preceding topol- 
ogy. In that case, a single Q value can be used to characterize an entire n-parton topology, irrespective 

'This distinction comes about from using quantities that are similar to a squared mass, squared transverse momentum, and 
squared energy, respectively. 
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Figure 3: The smooth-ordering factor (solid) compared to a strong-ordering function (dashed). 



of which branching was the last to occur. Since the distinction between Markovian and non-Markovian 
shower restart conditions only enters starting from the 4 — 5 parton evolution step, it will not affect our 
discussion of the second-order 2 — 3 branching process. For completeness, we note that the strongly 
ordered showers in VINCIA are of the ordinary non-Markovian type, while the smoothly ordered ones 
are Markovian. 



2.4 Smooth Ordering 

In addition to traditional (strongly ordered) showers, we shall also consider so-called smooth ordering 
[18]: applying the ordering criterion as a smooth dampening factor instead of as a step function. This 
is not as radical as it may seem at first. Applying a jet algorithm to any set of events will in general 
result in some small fraction of unordered clustering sequences. This is true even if the events were 
generated by a strongly ordered shower algorithm, if the jet clustering measure is not strictly identical to 
the shower ordering variable. An example of this, for strong ordering in p± and in m^i, clustered with 
the kr algorithm, can be found in [47]. 

In smooth ordering, the only quantity which must still be strictly ordered are the antenna invariant 
masses, a condition which follows from the nested antenna phase spaces and momentum conservation. 
Within each individual antenna, and between competing ones, the measure of evolution time is still 
provided by the ordering variable, which we therefore typically refer to as the "evolution variable" in 
this context (rather than the "ordering variable"), in order to prevent confusion with the strong-ordering 
case. The evolution variable can in principle still be chosen to be any of the possibilities given above, 
though we shall typically use 2p± for gluon emission and rrigg for gluon splitting. 

In terms of an arbitrary evolution variable, Q, the smooth-ordering factor is [47] 

^'"^^ ^ 0^ ' ^^^^ 

where Q is the evolution scale associated with the current branching, and Q measures the scale of the 
parton configuration before branching. A comparison to the strong-ordering step function is given in 
fig. 3, on a log-log scale. 
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In the strongly-ordered region of phase-space, Q <^ Q, v/e may rewrite the Pimp factor as 

Pimp - i + g - 1 Q2+---- (32) 

Applying this to the 2 — 3 antenna function whose leading singularity, eq. (16), is proportional to l/Q^, 
we see that the overall correction arising from the Q'^/Q'^ and higher terms is finite and of order l/Q"^; 
a power correction. The LL singular behaviour is therefore not affected. However, at the multiple- 
emission level, the l/Q^ terms do modify the subleading logarithmic structure, starting from 0{a1), as 
we shall return to below. 

In the unordered region of phase-space, Q > Q,\n& rewrite the Pimp factor as 




Pimp — ^9 Ao — ^9 I ^ ^9 ~l" • • • 1 ) (33) 



imp 

which thus effectively modifies the leading singularity of the LL 2 — )■ 3 function from l/Q^ to 1/Q^, 
removing it from the LL counting. The only effective terms oc \/Q^ arise from finite terms in the 
radiation functions, if any such are present, multiplied by the Pimp factor. Only a matching to the full 
tree-level 2 — )• 4 functions would enable a precise control over these terms. Up to any given fixed 
order, this can effectively be achieved by matching to tree-level matrix elements, as will be discussed in 
section 3.1. Matching beyond the fixed-order level is beyond the scope of this paper. We thus restrict 
ourselves to the observation that, at the LL level, smooth ordering is equivalent to strong ordering, with 
differences only appearing at the subleading level. 

The effective 2 — 4 probability in the shower arises from a sum over two different 2 — 3 (g) 2 — )• 3 
histories, each of which has the tree-level form 



A Pimp A OC J- ^-^ = , (34) 



thus we may also perceive the combined effect of the modification as the addition of a mass term in 
the denominator of the propagator factor of the previous splitting. In the strongly ordered region, this 
correction is negligible, whereas in the unordered region, there is a large suppression from the necessity 
of the propagator in the previous topology having to be very off-shell, which is reflected by the Pimp 
factor. Using the expansion for the unordered region, eq. (33), we also see that the effective 2 — )• 4 
radiation function, obtained from iterated 2 — )• 3 splittings, is modified as follows. 



in the unordered region. That is, the intermediate low scale Q, is removed from the effective 2 — >• 4 
function, by the application of the Pimp factor. 

Finally, to illustrate what happens to the Sudakov factors, we illustrate the path through phase space 
taken by an unordered shower history in fig. 4. An antenna starts showering at a scale equal to its 
invariant mass, ^/s, and a first 2 — ;> 3 branching occurs at the evolution scale Q. This produces the 
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Figure 4: Illustration of scales and Sudakov factors involved in an unordered sequence of two 2 — )• 3 
branchings, representing the smoothly ordered shower's approximation to a hard 2 — )• 4 process. 

overall Sudakov factor A2->3(-v/s, Q). A daughter antenna, produced by the branching, then starts 
showering at a scale equal to its own invariant mass, labeled ^/sl. However, for all scales larger than 
Q, the Pimp factor suppresses the evolution in this new dipole so that no leading logs are generated. To 
leading approximation, the effective Sudakov factor for the combined 2 — )• 4 splitting is therefore given 
by 

A1^4 ~ A2^3(^/^,Q) , (36) 

in the unordered region. Thus, we see that a dependence on the intermediate scale Q still remains in the 
effective Sudakov factor generated by the smooth-ordering procedure. Since Q < Q in the unordered 
region, the effective Sudakov suppression of these points might be "too strong". This situation is of 
course better than in the strongly ordered case, in which these points are zero, but it does suggest that 
a correction to the Sudakov factor may be desirable, in the unordered region. A study of Z — 4 jets at 
one loop would be required to shed further light on this question. 

Having presented introductory aspects of (antenna) showers, we now turn to a detailed discussion of 
how we match them to higher fixed-order calculations. 

3 Matched Antenna Showers 
3.1 Tree-Level Matching 

The strategy for matching to tree-level matrix elements used in VINCI A was defined by GKS in [18] and 
is tightly related to the veto algorithm outlined above. The philosophy is to view the shower produced 
by the smoothly ordered antenna functions as generating an all-orders approximation to QCD, which 
can be systematically improved, order by order, by including one more factor in the accept probability, 
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called the matrix-element correction. For a given trial branching, the full trial accept probability, up 
to the highest matched number of partons, is then obtained as a product of the ordinary trial-accept 
probability in the shower, multiplied by this extra correction factor. 

Since the shower is already correct in the soft and collinear limits, the matrix-element correction 
factor will tend to unity in those limits, but it can deviate on either side of unity outside those limits. As 
long as the combined accept probability is still smaller than unity, a probabilistic accept/reject step can 
still be applied, without having to worry about reweighting the events (which would be required if the 
total accept probability should exceed unity). It is also important to define the factor only in terms of 
physical cross sections (here represented by LO matrix elements), which guarantees positivity. (Again, 
if it were allowed to become negative, one would have to introduce negative-weight events, but this is 
avoided in the GKS strategy as defined in [18]). 

As we have seen, the shower furnishes an all-orders approximation to QCD. The aim is, for each 
resolved parton/jet multiplicity, to match the answer provided by the shower to an, ideally, all-orders 
exact expression, by applying a multiplicative correction factor, schematically [18,32] 

Matched = Approximate ^^^^^ . (37) 

Approximate 

At tree level, we in fact know only the first term in the expansion of the numerator, and we therefore 
expand the shower approximation to the same level. For n partons, assuming the approximation has 
already (recursively) been matched to the preceding order. 

Exact |M„|2 (38) 

Approximate y^ffrj Ctj Atj -Paccept \Mn-i^f = T^g^ -Pimp ^Ari Aj |Mn_iJ^ , (39) 
j j 

where subscript "T" indicates trial quantities (cf. section 2.2), we have suppressed the dependence on 
phase-space points, <I>, and the subscript j in the (n — l)-parton matrix element indicates the configura- 
tion obtained by performing the inverse shower step that removes parton j from the n-parton state. 

The factors in eq. (39) are easy to calculate if a tree-level matrix-element (ME) generator is available 
to provide the |Mp factors. The total ME-corrected accept probability is thus simply eq. (19), 

LO |Af„P 
-^accept = -Paccept Rn = ^accept o d p ^ \T7 12 ' ^'^^^ 

/ , j iJj L-j J"^imp Ari |-'Wn-lj | 

As mentioned above, this factor should be positive and smaller than unity, in order to avoid having to 
reweight any events. In practice, we have found the trial function defined in eq. (16) to guarantee this 
for all processes we have so far considered, mainly consisting of Z — )• n and H ^ n partons. As shown 
in [18], it is also possible to absorb subleading-colour corrections into this matching factor in a positive- 
definite way, but since subleading colour goes beyond the scope of our study we do not reproduce the 
arguments here. 

The fact that these factors change the distribution of the final set of generated events to the correct 
matrix-element answer can be explained by following the steps of the algorithm and summing over all 
possible branching histories. We start from Born-level matrix-element events, and generate trial shower 
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branchings, for a trial approximation to the {B + l)-parton matrix element of: 



\M^;^\' = E9T^CTaT^\Mh?\\ (41) 

i 

with i running over all possible single -parton clusterings that correspond to allowed shower branchings. 
Applying the LO accept probability, eq. (40), changes this to 



9Ti Cn An -Pac° pt l^a 



LO|2 



= |mLOJ2_ (42) 

That is, summed over shower histories, numerators and denominators are designed to cancel exactly, 
leaving only the LO matrix element for B+1 partons, as desired. Due to the full phase-space coverage 
and explicitly Markovian nature of the smoothly ordered shower algorithm, this procedure is straight- 
forward to iterate for Born -i- 2, 3, etc partons^. 

To provide a connection with antenna subtraction, which will be useful when we extend to NLO 
matching below, we can rewrite the ratio in eq. (37) by a trivial rearrangement, 

^, , , , . Exact — Approximate \ 

Matched = Approximate 1 + — . (43) 

\ Approximate / 

The numerator in this equation is very similar to a standard antenna-subtracted matrix element, though 
we emphasize that our antennae are of course modified by the presence of the Pimp and PAri factors. 

Let us finally re-emphasize that since we apply the correction factor to the antenna functions them- 
selves, thereby modifying the probability for a branching to occur, the probability for a branching not to 
occur is also modified. These corrections will therefore also be present in the Sudakov factors generated 
by the corrected shower evolution. The fact that the correction factor, ii„, is unity in all LL singular 
hmits (since the shower functions are guaranteed to match the matrix-element singularities there) im- 
plies that the LL properties of the Sudakov factors are not affected by this modification. However, the 
tower of subleading logarithms is changed, for better or worse. While it is known that finite terms do 
not exponentiate our corrections here also include a more subtle aspect, namely a resummation of the 
subleading logarithms present in the higher-order matrix elements. At this level, however, we cannot be 
sure that this procedure produces the correct subleading logarithms of a formally higher-order resum- 
mation. Therefore, we view it at present merely as an interesting, and hopefully beneficial, side-effect 
of unitarity-based matching. The examination of formally subleading terms carried out in this paper is 
a first step towards a more rigorous study of these aspects. 



'That is not the case for ordinary strongly ordered frameworks, due to the presence of dead zones in phase space and to 
the generally non-Markovian shower restart conditions. For such algorithms, addition of events in the dead zones [48], with 
CKKW-like Sudakov-factor prescriptions for multi-leg matching [49,50], would presumably be necessary. 
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3.2 One-Loop Matching at the Born Level 



For the Bom level, at NLO, the GKS matching strategy [18] reduces to the POWHEG one [7-9]. We 
nonetheless begin by recapitulating the steps used to derive the one-loop correction to the Bom-level 
matrix element, in our notation. We then continue to higher multiplicities. 

As our basis for one-loop matching we take the tree-level strategy described in section 3.1. Since 
the corrections are applied as modifications to the branching probabilities, without adding, subtracting, 
or reweighting events, the total inclusive rate after tree-level matching to any number of partons, is still 
just the leading-order. Born-level one. By the same token, after one-loop matching, at the integrated 
level, the total NLO correction to the inclusive rate must therefore just be the NLO "K-factor", 

^NLO 

^NLO ^ ^ (44) 
'^inc 

For processes like Z decay, where the NLO correction has no dependence on the Bom-level kinematics, 
this is trivial to implement as an overall reweighting factor on the Born-level events, 

vr 

where we have introduced the notation V for the NLO correction term, anticipating a similar notation 
for the multileg case below. Note that one could equally trivially normalize to NNLO or to data, as 
desired for the application at hand (we note though that such a normalization choice does not, by itself, 
ensure NNLO precision for any quantity besides the total inclusive rate). 

However, for more complicated cases, where the NLO correction does depend on the Bom-level 
kinematics, it is worth illustrating the general procedure for deriving a fully differential A'-factor, that 
can be applied to the Born-level events, for each phase-space point. This also serves as a useful warm-up 
exercise for the multi-leg case below. 

At NLO, we may distinguish between inclusive and exclusive rates for the first time. Either can in 
principle be used to derive matching equations between showers and fixed-order calculations, but the 
exclusive one is best suited for deriving expressions at the fully differential level. We recall that the 
exclusive n-jet cross section is defined as the cross section for observing n and only n jets, while the 
inclusive n-jet cross section counts the number of events with n or more jets. One therefore has the 
trivial relation 

^Z\Q) = Y.^T\Q)- (46) 

fc>n 

with Q the resolution scale of whatever (IR safe) algorithm is used to define the jets. 



3.2.1 Inclusive Born 

The total inclusive rate produced by the tree-level matched shower is just the Born-level matrix element, 

Approximate jM^ |^ , (47) 

where the subscript indicates the parton multiplicity (2 for Z ^ qq decay) and the superscript indicates 
the loop order beyond the Bom level (0 indicates the Bom loop order). Because cancellation of real and 
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virtual corrections is exact in both the unmatched shower as well as in the tree-level matching scheme 
described above, there are no further corrections to consider for the inclusive rate. In other words, the 
total integrated cross section produced by the shower is obtained merely by integrating eq. (47) over all 
of the Bom-level phase space. We now seek a correction term, V2, such that 

Matched ^ {1 + Vzz) \M^\'^ (48) 

gives the correct inclusive NLO rate. From eq. (45), we know that the correction term for Z decay is 

V2Z = - . (49) 

vr 

A systematic way of deriving this result, which can be applied to arbitrary processes, is provided by 
considering the cross section at the exclusive level. 



3.2.2 Exclusive Born 

The shower expression for the exclusive Z ^ qq rate (defined at the hadronization cut-off, which is the 
lowest meaningful resolution scale in the perturbative shower) is 

|m0|2 A(s,QLd) = l^2T (1-/' d<P.nt 9! CAg/gg + 0{al)] , (50) 

where we have expanded the Sudakov factor A to first order. Due to the presence of the hadronization 
scale, this expression is IR finite and can be defined in 4 dimensions. The colour factor for qq — )• qgq is 

and we assume that the antenna function. A, is either the one derived from Z decay [14] or has been 
matched to it, using LO matching. That is. 



\Mll 



g^, 2Cf Ag/^-^ = . (52) 



We first consider the limit Qhad — ^ 0, in which case the expression becomes 

IM2T A{s, 0) = |mO|' (^1 - d$ant g^s '^Cf Ag,^^ + 0{al)^ , (53) 

which can only be defined in the presence of an IR regularization scheme. We shall here use dimensional 
regularization, working in d = 4 — 2e dimensions. Below, we rederive the matching equations in 4 
dimensions, for Qhad / 0, and show that the same final matching factors are obtained in both cases. 
At NLO, the exclusive Z ^ qq rate at "infinite" perturbative resolution is 



\M^\' + 2Re[M'^M^ ] = \M^\' ( 1 + ^j^^U^ ) ' (54) 



19 



where we have written the right-hand side in a form similar to eq. (53), in d dimensions. Because the 
resolution scale has been taken to zero, there are no unresolved 3-parton configurations to include. The 
virtual matrix element is 



|M0|2 27r 



2Cf (2/,,-(e,/iV"i|)-4) , (55) 



with the function Iqq used to classify the e divergences [16, 27, 51]. Note that we have modified the 
definition of / to make it explicitly dimensionless, see appendix A. On the shower side, the integral of 
the Z — )• qgq antenna in eq. (53) is [16] 

^d^ant 2Cf gl Ag/qq = g 2Cf (-2/g5(e, /xVm|) + ^) , (56) 

and, not surprisingly, the difference comes out to be exactly ag/i^ x 1-^2 Writing this correction as a 
multiplicative K-factor, we obtain eq. (45). 

As a cross-check, we now repeat the derivation in 4 dimensions, reinstating the hadronization scale. 
The fixed-order side is then 

' I Jo ^''^^ I ' 

where the integral that has been added corresponds to unresolved 3-parton configurations, with A again 
given by eq. (52). Though eq. (50) is now defined entirely in 4 dimensions, we still need dimensional 
regularization to regulate the two last terms in the fixed-order expression. In principle, the integral in 
the last term could be carried out explicitly, but it is simpler to rewrite it as 

/ ' d«>ant 9s C Ag/gg = / d^ailt gS C Ag/gg " / d^ailt 5s ^ Ag/gg (58) 

^0 Jo JqI__, 



where the first term is just the full antenna integral, eq. (56), and the second term is identical to the 
one appearing in eq. (50), with which it cancels completely, cf. the definition of the tree-level matching, 
eq. (52). The final correction term is therefore again exactly equal to Os/tt x p. 

Note that the scale and scheme dependence of the Os/tt correction is not specified since its ambiguity 
is formally of order a^. For definiteness we take the renormalization scale for this correction to be 
proportional to the invariant mass of the system, fin = k™'^\/^ (so that /i^ = k^^'^mz at the Z pole), 
with A;™'^ thus representing the free parameter that governs the choice of renormalization scale for the 
total inclusive rate for Z — hadrons. We shall consider both one-loop and two-loop running options. 
The value of agimz) will be determined from LEP data in section 5. 

3.3 One-Loop Matching for Born + 1 Parton 

The approximation to the 3-parton exclusive rate produced by a shower matched to (at least) NLO for 
the 2-parton inclusive rate and to LO for the 3-parton one, is 

Approximate ^ (1 + V2) \mI\^ A2(m|, Qj) A^iQ^^, Qhad) , (59) 
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Figure 5: Illustration of the evolution scales and Sudakov factors appearing in the exclusive 3-jet cross 
section, eq. (59). 



where Mg is the tree-level Z — )• qgq matrix element and Qrz denotes the "restart scale". For strong 
ordering, Qrz is equal to Qs, while, for smooth ordering, it is given by the nested antenna phase spaces, 
i.e., by the successive antenna invariant masses. The subscripts on the two Sudakov factors A2 and A3 
make it explicit that they refer to the event as a whole, see the illustration in fig. 5. Again, we have the 
choice whether we wish to work in 4 dimensions, with a non-zero hadronization scale, Qhad> or in d 
dimensions with the hadronization scale taken to zero. We have maintained the hadronization scale in 
eq. (59), though we shall see below that the dependence on it does indeed cancel in the final result. 
The 2-parton Sudakov factor, A2, is generated by the (matched) evolution from 2 to 3 partons, 

2 

A2(m|, Ql) = \ - d$ant gl 2Cf Ag,gg + 0{al) , (60) 
JqI 

with Agjqq again defined by eq. (52). Notice that the integral only runs from the starting scale, m\, 
to the 3-parton resolution scale, Q\, hence this integral is IR finite, though it does contain logarithms. 
In the remainder of this paper, we shall work only with the leading-colour part of the Sudakov and 
matrix-element expressions, hence from now on we replace 2Cf in the above expression by Ca, 

A^'^Cml, Ql) = 1- I d$ant g's Ca + 0{al) . (61) 

The 3-parton Sudakov factor, A3, imposes exclusivity and is given by 

A3(m|, gLd) = 1 - E f[ d*^-* al {Ca Aej + 2Tr As,) + 0{al) , (62) 

where the index j runs over the qg and gq antennae, and we use subscripts E and S for gluon emission 
and gluon splitting, respectively. We have implicitly assumed smooth ordering here, which implies that 
the upper boundaries on the integrals are given by the respective dipole invariant masses (squared), sj. 
Note also that we must take into account all modifications that are applied to the LL antenna functions, 
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including Pimp, ^Ari, and LO matrix-element matching factors. (We do not write out these factors here, 
to avoid clutter.) I.e., the antenna functions in the above expression must be the ones actually generated 
by the shower algorithm, including the effect of any modifications imposed by vetos. 

For strong ordering, there are no Pimp factors, and the upper integral boundary is instead mm{Q'^, sj), 



2 

E 



min(Q|,Sj) 



d^ant g^s (Ca Aej + 2Tr Asj) + 0{al) 



(63) 



However, since strong ordering is not able to fill the entire 4-parton phase space [18, 23], full NLO 
matching can only be obtained for the smoothly ordered variant. It is nonetheless interesting to examine 
both types of shower algorithms, since even in the strongly ordered case, we may compare the Sudakov 
logarithms arising at 0{a1) to those present in the fixed-order calculation. 

On the fixed-order side, the expression for the 3-parton exclusive rate is simply 



Exact 



\M. 



0|2 



+ 



2Re[M3°M|*] 



+ 



Qhad d$4 



0|2 



(64) 



where the last term represents 4-parton configurations in which a single parton is unresolved with respect 
to the hadronization scale. For Z decay, d-dimensional expressions for the virtual matrix element have 
been available since long [16,26,27,52]. Details on the calculation and in particular its renormalization, 
are given in appendix B, in a notation convenient for our purposes. 

We now seek a fully differential matching factor, K3 = 1 + V3, such that the expansion of 



Matched = (1 + V3) Approximate , 



(65) 



reproduces the exact expression, eq. (64), to one-loop order. ("Approximate" here stands for the tree- 
level matched shower approximation, eq. (59).) Trivial algebra yields 
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JQI Jo 



QLd d$4 |M0|^ 



d$3 I Mil 2 



/ d^^ntgl {CaAej + 2TRAsj) 



(66) 



where we have reinstated d-dimensional forms of the one-loop matrix element and of the divergent 
3 — )• 4 terms. For a shower matched to IM4 ^ at leading order, the last two terms will cancel, by 
definition of the matched antenna functions (even for an unmatched shower, the difference could at most 
be a finite power correction in the hadronization scale, since the matrix element and the shower antenna 
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functions have the same singularities), yielding: 



T^LC 



LC 



V2Z 



f'TTl^ 2 p S ' 

+ / ' d$ant gs '^Ca Ag/gg + ^ / ^ cl ^ 9s {Ca Aej + 2Tr Agj) . (67) 
JQI Jo 

Rewriting the remaining integrals in terms of a set of standardized antenna subtraction terms, writing 
out the ordering functions for gluon emission and gluon splitting, Oe and O5, explicitly, and denoting 
the ARIADNE factor for gluon splitting by P4, we arrive at the following master equation for the second- 
order correction to the 3 -jet rate: 
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(69) 



with the standardized Gehrmann-Gehrmann-de Ridder-Glover (GGG) subtraction terms defined by [16]: 

19 
T 
17 

y 

1 
' 2 

whose IR limits and integrated pole structures were examined thoroughly in [16,27,51], though we have 
here rewritten the IR singularity operators /^^^ in explicitly dimensionless forms, see appendix A. (The 
alphabetical labeling in eqs. (69) corresponds to the notation used in [16].) 

The first line in eq. (68) represents a standard antenna-subtracted one-loop matrix element, normal- 
ized to the Born level, with the standardized subtraction terms tabulated in eq. (69), and the additional 
finite term V2Z originating from the NLO matching at the preceding order; see section 3.2, eq. (49). 

The subsequent lines express the difference between the simple fixed-order subtraction carried out in 
the first line and the actual terms that are generated by a matched Markovian antenna shower. Physically, 
these terms represent the difference between the evolution of a single dipole (the original qq system) and 
evolution of two dipoles (the post-branching qgq system), with a transition occuring at the branching 
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scale Q3. As mentioned above, the Oej and Osj factors in the third line represent the ordering criterion 
imposed in the evolution, either strong or smooth. For smooth ordering, they are 



I-Oej = 1 



Ql 



(70) 



1-Osj = 1 



Q'ej + Q3 



,2 ' 



(71) 



^Iq + Ql 



,2 ' 



with Qej the evolution variable used for gluon emissions, while for strong ordering, the factor (1 — Oj) 
can be removed if the integral boundaries are replaced by [Ql, sj] (note: this replacement should only 
be done in the third line). 

The last term in eq. (68) is an artifact of the ARIADNE factor, PAri> which was introduced in sec- 
tion 2.2 and is applied to gluon-splitting antennae in VINCIA. Summed over the two "sides" of the 
splitting gluon, this produces the same collinear singularities as the standard gluon-splitting antenna, 
but in highly asymmetric configurations in which the splitting gluon is near-collinear to a neighbouring 
colour line, the ARIADNE factor produces a strong suppression, which improves the agreement with the 
tree-level 4-parton matrix element [24], and which then generates an additional logarithm. 

Notice that all but the 6 A terms are defined in terms of standarized antenna functions, and the 
corresponding integrals can be carried out analytically, once and for all. We give explicit forms for each 
of these terms, for each choice of evolution variable, in the following section. 

The only terms of eq. (68) that need to be integrated numerically are thus the 6 A terms, which 
express the difference between the standardized antenna functions and those generated by the actual 
(matched) shower evolution, which may have different finite terms and/or be matched to the LO 4- 
parton matrix element. Nonetheless, since the previous lines already contain most of the structure, we 
expect these functions to be relatively well-behaved and numerically sub-leading. Specifically, the 5 A 
terms for gluon emission and gluon splitting, respectively, are defined by 



with A^^ the unmatched shower antenna function (as defined in [24, 33]) and the second-order LO 
matching factors, R^e and R^s (for ^ — ^ Q99Q and Z — )• qq'q'q, respectively), defined as in eq. (40), 
but including only the leading-colour terms in i?4g. For strong ordering, similarly to above, the Oj 
factors can be removed by changing the integration boundaries of the 5 A terms to [0, Q^]. 

Finally, we note that one could in principle equally well have defined eq. (68) without the terms on 
the third Une. The 5 A terms in eqs. (72) and (73) would then Ukewise have to be defined without Pimp 
and Paii factors. However, while this would give a seemingly cleaner relation with standard fixed-order 
subtraction, the behaviour of the (numerical) integrals over the 5 A terms would be more difficult, due 
to over-subtraction in the unordered regions. (Showers without either a strong-ordering condition or a 
smooth-ordering suppression greatly overestimate the real-radiation matrix elements in the unordered 
region [18,22,24].) Numerically, it is therefore more convenient to integrate the contributions repre- 
sented by the third line in eq. (68) analytically, leaving only the suppressed terms in eq. (73) to be 
integrated over numerically. 




(72) 



(73) 
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To be specific, the numerical integration over the SA terms is performed by rewriting the SA integrals 
in dimensionless MC form, as: 

2 N 

j=i i=i 

and similarly for the gluon-splitting terms, with $j a number of random (uniformly distributed) antenna 
phase-space points. The common factor 1/4 arises from combining the prefactor Svr^ above with the area 
of the phase-space triangle, 1/2, and the factor 1 / (167r^) from the phase-space factorization, d<I>ant- For 
smooth -ordering with an arbitrary normalization factor N± (so = N±p\), the ordering factors, 
Oj, reduce to: 

OEiqiQjAk^ qaghQcQk) = , (75) 

OE{qi, gjQk Qi, 9agb(ic) = same with z -H- A; , (76) 

Osiqi9j,qk ^ QaQ'bqc,q) = ^-^Vjk — ^ 

N±yjk + Xbc 

Osiqi, gjqk qi, qcq'tqa) = same with k , (78) 

where we have used y with ijk indices for the scaled invariants in the original qgq topology and x 
with abc indices for the integration variables in the antenna phase space. Note also that the y values 
are normalized to the full 3-parton CM energy (squared), while the x values are normalized to their 
respective dipole CM energies (squared). 



3.4 The Renormalization Term 



A further ingredient to be discussed is the choice of renormalization scale on both the fixed order and 
parton shower sides of the calculation, as these scales are in general chosen differently in both sides. 
Hence a translation term arises at second order, which must account for this difference, keeping in 
mind that, as the scale evolves from one to the other, flavour thresholds are passed. Our aim is to have 
the flexibility to use fixed order matrix elements renormalized according to their usual scheme, while 
maintaining the freedom to make a different choice on the shower side. 

The fixed order calculations for Z-decay to jets to which we match are customarily renormalized in 
a version of MS scheme called the Zero-Mass Variable Flavour Number Scheme (ZM-VFNS). In this 
scheme the bare QCD coupling is renormalized as 



9b = iJ-^giiA) 



1 + 



a 



Stt 



--+7E-ln47r + ln^ ) Pq 



(79) 



with ^0 = (IICa - 2nF)/3 = and ni? is the number of light flavours. One thus ignores flavours 
that are heavier than the scale of the calculation, both in the virtual and in the real calculations. Once all 
the UV poles are cancelled, one has a running coupling that depends on the number of light flavours for 
the scale at hand. One then changes the flavour number when a threshold is crossed. For our present 
case of Z boson decay to jets we take np = ^for jiR not too different from the ^'-boson mass. 
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Let us be more specific about the matching of as across flavour thresholds. At one loop, 



The value of Ap depends on the number of active flavours, as follows. When passing flavour thresholds 
the following one-loop matching conditions are imposed 

af\mb) = a(4)(mb), a^^\m,) = af\m,) . (81) 

These conditions can be satisfied if A^? obeys the matching conditions 

In = ^ In -rf^ ■ (82) 

With these conditions one can also express values for different flavour numbers into eachother. E.g. 

if rric < fJ-R < mi,, one can express ai'^\iJ-B) in terms of ai^'*(/iij) by the relation 

4'\^^r) = a^J^Hf^R)—-. (3), , • (83) 

For completeness we briefly review how this nj? -dependent UV singularity occurs in the context of 
the (inclusive) 3-jet rate, in the case where we only consider massless quarks [26, 52]. In the virtual 
contribution, the only one-loop diagram for Z — )• qqg that is sensitive to the number of flavours is the 
quark self-energy correction on the external gluon. The self-energy diagram itself, being scaleless, is 
zero in dimensional regularization. However, renormalization of the coupling amounts to adding a nj? 
counterterm on the exteral gluon line proportional to 

The real contribution contributes a nj? dependent (coUinear) 1/e pole as well, from gluon-splitting 

-„4(^)\ (85) 

In the sum over real and virtual contributions the poles cancel, as guaranteed by the KLN theorem, 
leaving a logarithm of the form 

On the shower side a related prescription is used, in which the running coupling is evaluated at a 
shower scale /ips, such that the scale again depends on the number of flavours. Depending on the value 
of /xps, a corresponding value of ni? is chosen, as well as of the QCD scale Kp- This is often different 
from that for a fixed order calculation. To give a specific example, matrix elements will typically be 
renormalized at a scale characteristic of the total CM energy, i.e., /i^p. = s an event-independent value, 
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while resummation arguments imply one best chooses a running scale, such as /ips = p±, for shower 
applications [28,29], which can differ per event. 

Shifting to a different scale for as of a given flavour number is quite straightforward. Translating 
from a shower scale /xps to a matrix-element scale /xme amounts to replacing, for an antenna function 

a,/ J + a,H^^^^^ In ( ^) + ©(a^)^ .1 . (87) 

For coherent parton-shower models, the arguments presented in [29] also motivate a change to 
a "Monte Carlo" scheme for Og, in which Aqcd is rescaled, for each np, by the so-called CMW 
factor ~ 1.5 (with some mild flavour dependence), relative to its MS value. If the shower model 
being matched employs this scheme, then a further rescaling of the renormalization-scale argument, 
Z^PS — ^ Aips/^CMW. should be used in eq. (87), with 



f67-3TT^-10np/3\ 
^^^^ = [ 2(33 - 2np) ) 



1.513 np = 6 

1.569 np = 5 

1.618 np = 4 ■ ^^^^ 

1.661 np = 3 



The translation of renormalization scale (and scheme) yields then an additional term to be added to the 
definition of V3 in eq. (68), 

= _<^]lNc_2np^J^\ ^ _^|i,(/4eA . (89) 
27r 6 V /^Ps / 27r 2 V /"ps / 

By inserting the above term, which enters at overall order ln(/Lf^g/ /Xpg), the two calculations can be 
compared consistently at one-loop accuracy. 

Note that if several different shower paths populate the same fixed-order phase-space point, then 
each path will in general be associated with a distinct /ips value. Thus, one Va^ term arises for each 
shower path, weighted by the relative contribution of each path to the total. Since for our case there is 
only one antenna contributing to Z — )• qgq, this particular complication does not arise here. 

We finally alert the reader regarding the use of different flavour number 's in the master equation 
(68). In that equation cancellation of 1/e divergences take place, already in the first line of the right 
hand side. For this cancellation it is important that the subtraction terms, originating from the shower 
expansion and listed in eq. (69), use as renormalized as in the fixed order calculation. All subsequent 
terms in the master equation are finite, and constitute differences of unordered and strongly ordered 
shower based terms, which are also finite, and beyond NLO. 
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3.5 Leading- Colour One-Loop Correction for Z — )■ 3 Jets 



Combining the results above, in particular eqs. (68), (69), and (89), we obtain the complete expression 
for the leading-colour^ one-loop correction for Z — 3 Jets, 
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(90) 



where: 

• the first line contains the full (leading-colour) one-loop matrix element, the V2Z correction from 
one-loop matching at the preceding order, and the Vs^ term from the choice of shower renormal- 
ization scale; 

• the second line contains the standardized subtraction term arising from the qg — qgg and gq — )■ 
ggq antennae; 

• the third line contains the standardized subtraction term arising from the qg — )• qq'q' and gq — )• 
q'q'q antennae; 

• the fourth to last lines contain the terms arising from the difference between the (matched) shower 
evolution and the standardized subtraction terms, including the consequences of ordering choices 
and modification factors such as those arising from the Ariadne factor and from matching to the 
LO matrix elements. 

We denote the singular subtracted 1-loop matrix element by S Virtual 



^We use the usual MC definition of leading colour and include terms oc Ca and ounp but neglect ones oc 1/ Ca 
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In section 4, we compute the analytical integrals corresponding to each of the shower-generated 
terms, for different choices of evolution variable, ordering criterion, and antenna functions. 

With the one-loop matrix element expressed as in appendix B.2, it is easy to see that the infrared 
singularity operators in eq. (91) cancel, leaving only explicitly finite remainders (which may still contain 
logarithms of resolved scales). This then constitutes the description of the one-loop matching for Z — 3 
jets, having already discussed the case for two jets. In the context of eq. (46) we have now corrected the 
first two terms on the rhs to NLO accuracy. 



3.6 One-Loop Correction for Born + 2 Partons 

To illustrate how the formalism presented here generalizes to higher multiplicities, we take the case of 
the NLO correction to Z — )• 4 partons. For simplicity, however, we continue to restrict our analysis of 
the correction factor to the leading-colour level. At NLO, the exclusive Z — > 4 partons rate at "infinite" 
perturbative resolution (similarly to above) is 

Exact ^ |M^P + 2Re[M4°M]*] . (92) 

Labeling the 4 partons by Z — )• z, j, k, £, there are two possible antenna-shower histories leading to 
each 4-parton configuration, with j and k the last emitted parton, respectively. Those two contributions 
both enter in the definition of the tree-level 4-parton matching factor, 

r. ^ |M4°(i,J,fc,^)P 

' Aj/iK\Mi{I,K,i)\'^ + Au/jL\Ml{i,J,L)\^' ^ ^ 

such that their sum reproduces the full 4-parton matrix element. Note that a separate such factor is 
applied to Z — )• qggq and Z — )• qq'q'q, and that we have suppressed colour and coupling factors here, 
for compactness (we ignore the small, non-singular extra interference terms for the special case where all 
four quarks have the same flavour). The antenna functions. A, are understood to include all such factors, 
as well as any Pimp and Pari factors appropriate to the branchings at hand. For a general n-parton matrix 
element, the denominator contains one term for each possible clustering. 

Labeling the IK — )• ijk history by A and the JL — )• jkl one by B, the sum over the two histories 
yields 

3 

P4A4(g4,0) ^f-^4|M3"pA2(m|,Q?)A3(Qf,g?) n(l + 0, (94) 

aeA,B m=2 

where it is understood that a is an index, not a power, and the last product factor takes into account the 
NLO matching at the preceding multiplicities. Expanding the Sudakov factors to first order and using 
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the definition of tlie NLO correction factor at tlie preceding multiplicity, eq. (67), tliis becomes 



k 



d^>anti?5 A 



4-s>5 



a|2 



aeA,B 



which we can rewrite as 
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(96) 



where we again emphasize that the antenna functions are understood to include all relevant coupling. 
Pimp, and Pari factors. The first term represents the new subtraction that the shower generates at 4 
partons, while the second represents part of the NLO correction to the preceding multiplicity. For one 
of the histories (the one followed by the "current" event), this correction has already been evaluated and 
can be reused. The contribution from the other history will have to be recomputed, however In general, 
there will be one subtraction to perform at the n-parton level, and there will be m ~ n — nBorn — 1 new 
subtractions that have to be done at the (n — l)-parton level, in addition to the one that was already done 
during the evolution of the current event. 

Clearly, there is an undesirable scaling behavior built into this, which will make NLO matching at 
many partons quite computing intensive. An alternative, which eliminates the sum over histories, is that 
of sector showers, see e.g., [24, 53]. Though this is not the main avenue pursued in this paper, we there- 
fore give some comments below on how a sector-based NLO matching scheme could be constructed. 



3.7 One-Loop Matching for Sector Showers 

The matching conditions derived above may also be applied to so-called sector showers [24,53], with a 
few relatively minor modifications. The expansion of the Sudakov factors generated by the LO matched 
shower will now contain integrals over ratios of matrix elements (which are the LO matched sector 
antenna functions), multiplied by sector vetos. The presence of the sector vetos makes analytical phase- 
space integration more difficult. 

However, since the sector approach merely represents a different way of decomposing the same 
singularities as the global one, we may effectively recycle the integrals carried out for the global case 
by adding and subtracting the terms produced by a smoothly ordered "standard" shower (i.e., using the 
GGG functions). The first four lines of eq. (68) then remain unchanged. The definition of the 6A terms 
in the last line, however, changes to 




for the terms arising from the qgq — )■ qggq Sudakov factor, with analogous ones arising from the gluon- 
splitting contributions. The step function, Qj, represents the sector veto applied to the sector antenna 
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functions. The sector antenna function, up to the tree-level matched orders, is just 



9iC = . (98) 

Since the individual sector and global antenna functions have different singularity structures (they 
are only guaranteed to have the same singularities at the summed level), the integrals in eq. (97) are 
divergent, and cannot be carried out numerically. In order to obtain numerically convergent integrals, 
we must divide up the contributions of the global terms onto each sector, and perform a set of correlated 
integrals in which the singularities explicitly cancel in the divergent limits, 

2/^ r\sct I A set Astd \ I i/f, ^2/^ osct Astd 



where each line now corresponds to the sum of contributions to a single sector, for which the difference 
between sector and global antennae is finite. The individual integrals are of course still divergent, but 
they can now be treated numerically by collecting the terms on each line under a single integral sign. 
Analytically, this is complicated since the two integrals on each line are not associated with the same 
kinematics map*^. Numerically, however, we may still ensure a point-by-point cancellation in the sin- 
gular limits by keeping the two integrals formally separate, but carrying them out simultaneously, in a 
correlated way, as follows. 

For each antenna, generate a random uniformly distributed phase-space point, 

[Vij, Vjk, </>] , (100) 

and evaluate the first term in eq. (99). (If the point is outside the relevant sector, this term is zero for the 
time being). If the pair (jk), say, corresponds to a sector shared with a neighboring antenna (as in our 
example), check whether the correlated phase- space point defined by 

[Vjk, 1 - Vjk - Vij, </>] , (101) 

in the neighboring antenna passes the same sector veto as before, and if so, subtract the global term 
corresponding to the second term on the first line of eq. (99). The sampling represented by eq. (101) is 
uniform, as long as that used to generate the original point, eq. (100), is uniform. The replacement of 
Uij by 1 — Ujk — Vij corresponds precisely to the mapping z — )• 1 — z in the collinear limit, which is 
what is required to reestablish a point-by-point cancellation of the sector and global singularities. 



4 Sudakov Integrals 

In this section, we work out the standardized Sudakov integrals appearing in the second and third line of 
eq. (68), for each choice of evolution variable. We also study the soft and collinear limits of the Sudakov 

'in the example here, for Z — > 3 partons, one term is associated with branchings in the qg antenna, while the other is 
associated with branchings in the gq one. 
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integrals and compare them to those of the one-loop matrix element. This provides an explicit check of 
whether the first-order expansion of the Sudakov factors generates the correct logarithms present in the 
fixed-order calculation. 

Given our choice of the GGG antenna functions as our standard ones, the relevant terms are 



9l 



S 9 9 

(102) 

The general form of the first term, which originates from the 2 — )• 3 branching step, is 

alCA j 4 d$ant = ^ (E ^^^^(^' ^i)) (103) 



where the definitions for the Ki and the /j functions are given in appendix C, for each type of antenna 
function and ordering variable. Their derivation and soft/coUinear structure will be discussed more 
closely below, for each choice of ordering and evolution variable. The form of the 3 — )• 4 integrals 
depends on whether we work in the context of strong or smooth ordering. We shall now consider each 
of those cases in turn, beginning with strong ordering. 



4.1 Strong Ordering 

For strong ordering, the inverted ordering conditions in eq. (68), (1 — Oe^/Sj)' reduce to step functions 
expressing integration over the unordered region. The integration surface is thus limited from below by 
the phase-space contour defined by the evolution scale of the first branching, Q^, and from above by the 
edge defined by the invariant mass of the antenna. 

The expression generated by the 3 — >• 4 splitting case for gluon emission is 

(104) 

where Ki and li are the same as those for the 2 — )• 3 term above, though they here appear with different 
arguments. The remaining case is the 3 — )• 4 gluon splitting defined by 

2 

-g^^Y^^FPA, r\l-Os,)eld^,^, = -^P^^^His,„Ql) -^P^^^His,„Ql). (105) 
j=l -^0 

with H defined in appendix C and Pa^ as defined in eq. (20). We will discuss the derivation of these 
terms in more detail in the following three subsections, for strong m/)-, p±, and energy-ordering, re- 
spectively. 



4.1.1 Dipole Virtuality 

We begin with dipole virtuality as evolution variable, which is perhaps the simplest case. We start by 
repeating the integrals of eq. (68) with the one-particle phase space defined as in eq. (8). In the case of 
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dipole virtuality the contours are triangular (fig. la). We recall that, for the g ^ qq terms, it is the qq 
invariant mass that is used as evolution variable, regardless of what choice is made for gluon emissions. 
The m£) scale of the previous emission still enters, however, since that defines the ordering scale applied 
to both emissions and splittings. The explicit forms of the terms in eq. (102) are: 
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(106) 

as defined in eq. (20) and the gluon-splitting antenna 63 has its 



with Pa 

3 cqg-rsgq 

singularities in si. 

For compactness, we only show the integration for the double-pole (soft-coUinear eikonal) terms 
present in both Og and dg here, which are the only sources of transcendentality-2 terms. The full antenna 
integrals, including also the lower-transcendentality terms originating from single poles and finite terms, 
are given in appendix C. The T = 2 part of the integral is 
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(107) 



lmm{Sqg,Sgq) 

To evaluate this expression, we first rewrite it in a dimensionless form in terms of the rescaled integration 
variables Ui = Si/{s — with upper boundary 1 and lower boundary 

mm{Sgg,Sgg) 



(108) 



mm{Sgg,Sgg) ' 

The integration is over a triangular surface. The lower integration boundary cuts off the evolution 
variable at the value of the 3-parton evolution scale. The other boundary is determined by the total 
energy of the dipole before branching which here is s. We use the integral 

dy. 
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4 Sudakov terms, let us for definiteness assume that we are in a 3-parton phase- 
The opposite case is symmetric. Again we only include the T = 2 terms 



To discuss the 3 - 
space point with Sqg 

exphcitly here, with the details of the full antenna integrals relegated to appendix C. 
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The integration is again over a triangular surface. The total energy of the dipole before branching is 
now Sqg. The integral in eq. eq. (1 1 1) corresponding to the sum over antenna integrals only contains one 
dg integral because the other has equal upper and lower integration boundaries. Note that this integral 
actually vanishes for Sqg < Q^, which amounts to the dipole- virtuality ordering allowing the 3 — >• 4 
branchings to populate their full respective phase spaces (i.e. no correction term is necessary). 

Focusing on the case Sqg > 2sgq for which the second integral is nonvanishing (which now amounts 
to the ordering condition imposing a nontrivial restriction on the 3 — > 4 phase space, see fig. 2a), we 
obtain, including the 2 — 3 term 
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with lower-transcendentality terms again available in appendix C. For the mirror case 

tially symmetric expressions are obtained, while for the intermediate cases in which the two invariants 

are within a factor 2 of each other, the second integral in eq. eq. (112) simply vanishes. 

The full double-logarithmic term from the expanded Sudakov terms in eq. (106), for strong ordering 
in dipole virtuality, is then 
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where the © function ensures that the second term is only active if 
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so that the expression applies over all of phase space. 

We shall now consider the infrared limits of this result, and compare them to those of the one-loop 
matrix element. For this comparison we keep only terms that involve logarithms of the invariants. The 
soft limit corresponds to vanishing Q\ (.^min — ^ 0). The first line of eq. (114) represents the contribution 
of the 2 — 3 expanded Sudakov. To find the contribution in the soft limit, we choose to approach the 
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limit along the diagonal of the phase space triangle. Parametrizing this by Sqg/ s 
for this term 



y we find 
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The contributions of the 3 — 4 Sudakovs in the soft limit are examined in two separate cases corre- 
sponding to the two regions in fig. 2a. In the first case given by Smax < 2smin, con^esponding to the 
light grey area in the figure, the step function in eq. (114) yields zero. In the second case given by 
•Smax > 2smin, Corresponding to the dark grey area in the figure, the step function is equal to one. The 
double logs and dilogarithms now yield a finite contribution that does not diverge in the soft limit. We 
can understand this by parametrizing the soft limit by A 
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so that the integral becomes 
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This implies that the integration variable scales with the integration limits and is independent of the soft 
limit. We can also expect this behaviour from examing fig. la. The shape of the different regions does 
not change for different values of Q^, in contrast with the case of transverse momentum, as we will see 
below. 

After the poles cancel in eq. (90), the pole-subtracted version of the one-loop matrix element, 
S Virtual , defined in eq. (91), contains all the relevant terms on the matrix-element side. The transcendentality- 
2 terms of SVirtual are given by 
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Including the transcendentality-1 terms (see appendix B), taking the soft limit by sending Sqg/s = 
Sgg/s = y — 0, and keeping only logarithmic terms, the pole-subtracted matrix element (ME) reduces 
to 
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The single log proportional to Ca originates from the renormalization term and the single log of the 
closed quark loops (proportional to n f) arises due to the definition of the infrared singularity operator, 
defined in the appendix in eq. (159). 

Taking the same limit of the Sudakov integrals for dipole virtuality eq. (106), but omitting for the 
time being the renormalization term, Vs^, we find for the parton shower (PS), 
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We see that the soft limit almost cancels against eq. (119). For an NLL-accurate shower, however, 
aU divergent terms should match precisely, leaving at most a finite remainder in the final matching 
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correction, eq. (90). In the expressions above, this holds for the In^(y) term but not for the single 
logarithms (different coefficient). Interestingly, the remainder is proportional to the QCD /3 function, 
specifically 
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It can therefore be absorbed in the choice of renormalization scale by solving for /ips in Va^, which 
yields: 



MPS « y s • 
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This tells us that, in the soft limit, the specific choice of a renormalization scale that is hnear in the 
branching invariants will absorb all logarithms up to and including ln(y). Interestingly, this reasoning 
would rule out fi'j^ (x p^, since our -definition is quadratic in the invariants, p'j_ = SijSjk/s. A better 
choice of renormalization scale would appear to be /zr oc m/j, specifically 
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Taken at face value, this seems to contradict the standard literature [28] on p± as the optimal universal 
renormalization- scale choice. However, as we shall see below in fig. 6, there is in fact no choice of 
renormalization scale that absorbs all logarithms for this particular evolution variable; the choice fiji oc 
m £, merely manages to reabsorb the additional logarithms that are generated by the ordering condition 
as y — 0, but leftover logs in other parts of phase space will remain uncanceled, ruining the NLL 
precision. In that sense, choosing fiji oc p± would simply leave a different set of uncanceled logs, 
nonvanishing as y — )• 0. 

Before we show the results over all of phase space however, we first investigate a complementary in- 
teresting limit, the hard-collinear one, which is characterized by one of the invariants becoming maximal 
while the other vanishes. In this limit, the pole-subtracted one-loop matrix element, SVirtual , becomes 
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There are no log-squared terms in this limit, and both of the single-log terms are half as large here as 
they were in the soft limit. 

The Sudakov integrals for m/) -ordering yield one divergent term, — \Ca ln(y), in the hard-collinear 
region, modulo a factor as/ (27r). The Sudakov integral for gluon splitting in the neighbouring antenna, 
represented by the first term on the second-to-last line of eq. (90) is specified in the last line of eq. (106). 
The step function is only non-zero for the first term in the hard-collinear limit Sqg — )• s, Sgq — )• and 
produces a term gPyi^ rii? ln(y). The numerator of the corresponding Ariadne factor contains the invari- 
ant of the neighboring dipole Sgq which vanishes in this limit and causes the dipole splitting contribution 
to reduce to zero. The nj? -dependent contribution is instead shifted to the last term of eq. (90), which 
has the same limit but without the Ariadne pre-factor. The hard-collinear limit of the shower terms, 
including only terms involving logarithms of the invariants and not including the Vs^ term, is therefore 
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(a) ^ps = x/i 
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(d) ^ps = imc, with CMW 
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Figure 6: NLO coirection factor for strong m/5-ordering, with GGG antennae. Top row: fiji = ^/s 
(left), = p± (middle), and /xr = nif) (right). Bottom row: using the CMW Amc> with fiR = \mD 
(left) and ij,r = mo (right). For all plots, Og = 0.12, np = 5, and gluon splittings were evolved in rugq. 



Again, the combination (ME — PS) relevant for computing the correction factor is proportional to the 
QCD /3 function, and in fact has exactly the same form as eq. (121). The conclusion is therefore that, 
also in this limit, all logarithms through ln(y) can be absorbed by choosing a renormalization scale 
which is linear in the vanishing invariant. The particular choice which is linear in both the soft and 
collinear limits is ^ps oc m/j. To illustrate this, we show the full NLO Z — 3 jets correction factors, 
(1 + Vsz), for m^i-ordering with a few different choices of renormalization scale and scheme, in fig. 6. 
Note that axes are logarithmic, in In(yy ) = ln{sij/s), to make the infrared limits clearly visible. 

Without the Vs^ term, the correction factor looks as depicted in the top left-hand plot in fig. 6. The 
increasing contours towards the axes indicate uncanceled logarithms in the correction factor. The middle 
pane shows the correction factor derived for /ips = p±- As discussed above, there is an uncanceled 
logarithm in the soft limit (lower left-hand corner of the plot), since p± is quadratic in the vanishing 
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invariants there. However, in the hard-collinear limits (upper left-hand and lower right-hand corners), 
p± is linear in the vanishing invariant, and hence the contours remain bounded there. In the right-hand 
pane, we show the choice ^ps = mo, which can be seen to lead to bounded correction factors (below 
~ 1.3) in all three phase-space comers. Nonetheless, there is still an uncanceled divergence between the 
soft and hard coUinear limits. We shall see in the section on -ordering below that the cure for this is 
basically to choose a better evolution variable. 

In the bottom row of fig. 6, we show a few variations on ^ps = rriD, specifically we include the 
CMW rescaling of Aqcd> as defined by eq. (88), and show how a variation of the renormalization scale 
by a factor of 2 affects the correction factor. In the left-hand pane, we show fips = ^rriD and in the 
right-hand one ^ps = m£). Of these, the choice /ips = ^m£,, with CMW rescaling, leads to the 
smallest correction factors (best LO behaviour), and this could therefore be taken as a useful default for 
mD-ordering, e.g. for uncertainty estimates. 



4.1.2 Transverse Momentum 

For a shower ordered in p±, the antenna phase-space integrals in eq. (68) are performed over contours 
such as those depicted in fig. le. The curved contours motivate a coordinate transformation from (si , S2) 
to a basis defined as the dimensionless evolution variable y = = , complemented by an energy- 
sharing variable, which we define as 2: = ^. Note that the coordinate transformation depends explicitly 
on the total invariant mass s of the 2 — )• 3 dipole. For the 3 — )• 4 integrations, the invariant mass s is 
replaced by the invariant mass of the appropriate dipole (either Sqg or Sgq). The integration boundaries in 
z are determined by the intersection of the invariant mass of the dipole with the evolution parameter Q^. 
The choice of y and its integration boundaries make the effect of strong ordering especially clear since 
we see integration from to the total invariant mass of the dipole (the unordered region). As before, 
the integration over the gluon-splitting antenna (eg) makes use of a different phase space integration, in 
mqg, and only uses the evolution parameter as a cut-off in the singularity of the corresponding dipole. 
The contributing terms are: 
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For n 



1 we set rrij^ 
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Sgq. The Ariadne factor . 



is defined in eq. (20). The max condition on the outer integration boundary of A2 and ^3 reflect that 
the correction term disappears if the generated Q\ is larger than the invariant mass of the dipole. As for 
m/)-ordering, we here work out the most divergent behavior explicitly, by focussing on the double log 
terms arising from the eikonal term in the antenna, and relegate the full form of the antenna integrals to 
appendix C. The double poles give rise to terms 
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which lead to the following generic transcendentality-2 integrals. 
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The double logarithm in the shower expansion is generated by a combination of the 2 
Sudakov integrals, with the respective pieces adding up to 
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We see that a partial cancellation arises between the first two terms (which come from the 2 — )• 3 
Sudakov expansion) and the last three (which come from the 3 — > 4 expansion). What remains is a log 
squared in both invariants In {sqg/ s) In (sgg/s). 

At the single-log level, the 3 — )• 4 terms give a numerically larger coefficient than the 2 — > 3 ones, 
leading to a single log remainder. The gluon-splitting term also reduces to a single log. The overall 
result in the soft limit is then 
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Comparing with the result of the virtual correction in the soft limit, eq. (119), we see that the shower 
generates the double log terms correctly, and, similarly to the case of m/j-ordering, there is a single- 
log remainder which is proportional to the QCD /3 function. However, for -ordering the constant of 
proportionality is 1, rather than ^, a difference which translates to the optimal renormalization-scale 
choice being quadratic in the invariants in this case, rather than linear. Before commenting further on 
this difference, let us first consider the complementary, hard-coUinear, limit. 

In the hard-coUinear limit, we find the same as for m^i -ordering. 
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Double logs (eikonal parts of the antenna) also appear at both the 2 — 3 and 3 — )• 4 levels, but cancel 
among each other as almost all other antenna terms do; what remains at the single-log level is the inte- 
grated difference between a quark-antiquark antenna and a quark-gluon antenna, plus the ni?-dependent 
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'Ariadne Log'. The only contributing Sudakov gluon splitting contribution is the second term in the last 
line of eq. (126). Integration over the Sgq dipole, however, is associated with an Ariadne factor carrying 
Sqg in the numerator and therefore reduces to zero. As before, we can write the remainder as half the 
QCD /3 function, which implies that a renormalization scale Unear in the vanishing invariants can absorb 
the logarithm. 

To summarize, for -ordering we find that the optimal renormalization-scale choice must be quadratic 
in the vanishing invariants in the soft limit and linear in the hard-collinear limit. Those conditions are 
fulfilled by p± itself, thus 

oc pi = ^ (133) 

absorbs all logarithmic terms up to and including a1 ln(y) in the LO couplings. 

Illustrations of the full NLO correction factors, (1 + V^z), are given in fig. 7. The ordering of 
the plots in the top row are the same as in fig. 7, showing, from left to right, /ips = -y/s, /ips = P^_, 
/ipg = m£,. Similarly to the case of strong m^) -ordering, both of the latter two choices exhibit no 
logarithmic divergences in the hard-collinear regions (top left and bottom right corners of the plots), 
but in the soft region (bottom left corner) it is here pps = p± which leaves the correction factor free 
of logarithms. Indeed, we see that the combination of evolution and renormalization in p± leads to a 
rather flat correction factor over all of phase space, showing that this combination is indeed "better" than 
m£)-ordering. 

In the bottom row of plots in fig. 7, we include the CMW factor and show the correction factors 
for fips = p± (left) and fips = 2p± (right). In particular on the left-hand pane, the NLO correction 
factor is essentially unity in the soft Umit, while the corrections in the hard-collinear regions remain less 
than ~ 20%. This gives some additional weight to the arguments for -ordered showers with p± as 
renormalization scale being the best default choice for strongly ordered dipole-antenna showers. It also 
provides some rationale why one typically finds a rather large value of as{mz) ~ 0.13 (with CMW 
rescaling, or as(m^) ~ 0.14 without it) when tuning such models to LEP event shapes; there is still a 
genuine order 20% NLO correction in the hard resolved region (upper right-hand corner). We return to 
this in more detail in the context of full LO -i- NLO matching in section 5. 



4.1.3 Energy 

To put the differences between mo and p± in context, we now briefly examine the case of energy 
ordering, which is known to produce the wrong DGLAP evolution in the coUinear limit [21,22,54], and 
hence we should find larger (possibly divergent) NLO corrections. 

Slicing phase space with the energy variable Q| = Sijk{yij + yju)"^, see fig. If, requires the use 
of an explicit infrared cut-off because the contours otherwise allow for the invariants to hit singular 
regions for every value of the contour. We will here use a cut-off in transverse momentum (a cut-off in 
dipole virtuality is also possible). The cutoff motivates us to switch to a different choice of integration 
variables. Therefore integration is transformed from (si, S2) to the dimensionless evolution parameters 
y|; = ^ = ^'^^'^^^ and completed with the energy sharing variable C, = . The interesting integrals 
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Figure 7: NLO correction factor for strong -ordering, with GGG antennae. Top row: ^ups = ^/s 
(left), ij,ps = p± (middle), and /ips = mn (right). Bottom row: using the CMW Amc. with fips = p± 
(left) and fips = 2p± (right). For all plots, as = 0.12, np = 5, and gluon splittings were evolved in 
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with 



'"/if 

With A£\ = 03, A£% = ^3 and A£^ = 63. The inner integral has been rescaled to make it independent 
of the outer integral with Q = ueC ■ To establish the cut-off, we use the relation = Ap^/s, which 
we demand to be larger than cut-off A. In terms of our variables we then have the condition 
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The upper and lower limits on C are then 
Focussing on the eikonal integral 
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the result for this integral is 
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Thus we see that there are explicit non-cancelling double-logarithmic terms that involve the hadroniza- 
tion cutoff, A. Depending on the ratio between the dipole mass and the cutoff, these would lead to 
asymptotically divergent correction factors. 

One might wonder whether using a linearized form of energy ordering would make a difference, 
see fig. Ic. Rather than go through the derivations again, we merely show the full NLO corrections in 
fig. 8, for both linear (top row) and squared (bottom row) energy ordering, for an (arbitrary) dimension- 
less cutoff value of A = 10~^. 

From left to right in both rows, we show the three renormalization-scale choices, ^ps = p± (left), 
/xps = mn (middle), and fips = Qe (right), with the latter equal to linear energy in the top row and 
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Figure 8: NLO correction factor for strong energy-ordering, with GGG antennae, for various 
renormalization-scale choices and linear (top row) and squared (bottom row) scaling of the evolution 
variable with gluon energy. 

squared energy in the bottom one. Interestingly, the contours in the linear case are increasing towards 
the soft region, while they ultimately decrease in the squared case. It is clear, however, that no intelligent 
choice of renormalization scale can absorb the infrared divergences. Moreover, any other choice of A 
would have led to different contours, due to the explicit ln(A) terms, hence even if a "least bad" choice 
was found, it would not be universal. 

As mentioned above, the main point of showing these comparisons is to place the comparison be- 
tween m D and p_L in the previous subsections in perspective. Thus, while we saw that -p^ was generating 
a better-behaved correction factor than m^^, the one for m^) is still far better behaved than is the case for 
energy ordering. From this perspective, we still believe it could make sense, e.g., to use m/j-ordering, 
with the NLO correction factor included, as a conservative uncertainty variation for a central prediction 
based on -ordering. 
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4.2 Smooth Ordering 



We will now discuss the same Sudakov integrals as in the previous subsections but for the case of smooth 
ordering (section 2.4). This is especially interesting given that smooth ordering is the way VINCIA is able 
to fill all of phase space without significant under- or overcounting at the LO level [18]. As discussed in 
section 2.4, however, this does involve some ambiguity in what Sudakov factors are associated with the 
unordered points, and the NLO 3-jet correction factors should tell us explicitly whether this ambiguity 
generates problems at this level. 

The Sudakov integrations are actually more straightforward for smooth ordering than was the case 
for strong ordering, since the Pimp factor regulates the integrands on the boundaries. Therefore the 
integrations always run over the full phase space of the system. The 2 — )• 3 splitting generates the same 
terms as in the strong-ordering case, eq. (103). Including also the 3 — )• 4 terms, the expanded Sudakov 
generates the following antenna integrals, 



9's 



Jo ~{ Jo ^Ej + ^3 ~l Jo rUgg + 



(140) 

where Qs is the evolution scale evaluated on the 3-parton configuration and Qe (rrigq) is the scale of 
the 3 — )• 4 emissions (splittings) being integrated over. The full answer for the 3 — )• 4 case for gluon 
emission is 



(141) 

where Ki and Li can be found in appendix C. The full answer for the 3 — )■ 4 case for gluon splitting is 



d E uePa, r el d$.„t = Ql) - ^Gis,„ QI). (U2) 

Jo lllgg i- ^TT ^TT 



"qq 



where G can be found in the appendix. We will discuss the derivation of these terms in more detail in 
the following two subsections, for smooth mo- and -ordering, respectively. 
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4.2.1 Dipole virtuality 



Since the 2 — > 3 emission terms remain the same as in the case of strong mD-ordering, we only need to 
rederive the 3 — >• 4 contributions to V^z, which are 







Qej + ^3 



Q 



47r 



+ - 



Ca 

Sqg \J0 

Ca 

^gq \Jo 



2 •'99 



dS2 



fqg-S2 



dsi 



2S2 



2 "31 



dS9 



S2 

Sgg-S2 



ds 



32 



Qi + 2s2 



+ 



2°ig 



dsi 



2si 



63 d$ant 



+ 



2 "99 



dsi 



2s 



/ /•■599 rSqg-S2 

+2nF — / ds2 / dsi + 



PA2 



Sgg~S2 



4 



4 



dsi 



Sl 



(143) 



with = 2 min(sqg, Sgq) and Cg carrying the singularity in si. We will focus again on deriving the 
transcendality-2 terms explicitly, with the full expressions given in the appendix. We start by recalling 
the expression for the strongly-ordered 2 — 3 branching, 
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where we have transformed yi 



for i = 1,2 and y| 
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soft region y| — 2 (since we take the invariants as vanishing simultaneously), we see that the remainder 
is just a finite term that contains no logarithms of the vanishing invariants. 
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We will receive this contribution twice. Including all divergent logarithmic contributions and disregard- 
ing constant terms such as in eq. ( 145) , we find the same as in the strong-ordering case. 
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Figure 9: NLO correction factor for smooth m/j -ordering, with GGG antennae, and /Ups = mo (left), 
^J'PS = ^f^D (middle), and /ips = ^thd with CMW rescaling (right). For all plots, = 0.12, np = 5, 
and the evolution scale for gluon splittings was niqq. 



and hence the preferred choice of scale in the soft limit remains one which is linear in the vanishing 
invariants, such as /ips oc mo- 

In the hard coUinear limit the Sudakov integrals plus the 'Ariadne Log' reduce to 
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again the same as in the strongly ordered case, cf. eq. (125). 

To summarize, we do not expect any qualitatively different limiting behaviour in the smoothly or- 
dered case with respect to the strongly ordered one, though details may of course still vary. To illustrate 
this, we include the plots in fig. 9. In all cases, we use a renormalization scale oc m/j, but with different 
prefactors, from left to right: /xps = m/j, /ips = mj:,/2, and finally ^ups = mr)/2 with CMW rescal- 
ing. In particular the latter gives correction factors very close to unity in both the soft and hard coUinear 
limits, while we still see the leftover divergence inbetween those limits that was also present in the case 
of strong m,£)-ordering, cf. fig. 6. Nonetheless, it is worth noting that for a large region of phase space, 
say with niij > 0.1m (corresponding to In(yjj) > —4.6), the corrections are still quite well behaved 
and relatively small, less than ~ 20%. 
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4.2.2 Transverse momentum 



Again we only need to recompute the contributions from the 3 — )■ 4 Sudakov terms, as the 2 — )• 3 ones 
are the same as in the strongly ordered case. The 3 — )■ 4 Sudakov integrals are 
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As before we focus on explicitly calculating the transcendentality-2 contribution arising from the eikonal 
part of the antenna in the first term in the first line of eq. (148), 
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where we have transformed yi 
y| — >• 0, this yields 
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Adding the contributions from the 2 — 3 splitting and transcendentality-1 terms, we find the following 
result for the soft limit 
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as in the strongly ordered case. The double logarithm matches with SVirtual and the single logarithm 
can be absorbed by choosing a renormalization scale that is quadratic in the vanishing invariants, such 
as ^ps oc p_L- 

In the hard coUinear limit, the shower integrals behave as 
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the same as in all the other cases. This completes the argument that indeed ^ps oc p_\_ is the appropriate 
choice also for smooth -ordering. 

In fig. 10, we show the NLO correction factors, (1 + V32), for smooth -ordering. The top row 
shows the correction factors without using the CMW rescaling of Aqcd> and the plots in the bottom 
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In(jij) ln(iij) 

(c) ^ps = p_L, with CMW (d) ^^ps = 2p±, with CMW 



Figure 10: NLO correction factor for smooth -ordering, with GGG antennae, without (top row) and 
with (bottom row) the CMW rescaling of Aqcd- The left-hand panes use /xps = p± and the right-hand 
ones /ips = 2p±. For all plots, ag = 0.12, np = 5, and the evolution scale for gluon splittings was 

mqq. 
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(a) ^ps = p±, N± = I (b) ij.ps = p±, N± = A, g ^ qq'm p± (c) ^ps ^ p±, N± =4, Paj = 1 

Figure 1 1 : NLO correction factor for smooth -ordering, with GGG antennae: variations of how gluon 
spUttings are interleaved with gluon emissions, see text. We used as = 0.12, np = 5, and /ips = p±. 

row include it. For the left-hand panes, we used a shower renormalization scale of /ips = p± , and for 
the right-hand ones we used /ips = 2p± . 

We see that all correction factors are essentially well-behaved, with no runaway logs, similarly to the 
case of strong -ordering. However, for the case of smooth -ordering, it looks as if the CMW rescal- 
ing (bottom row) is almost doing "too much" in the soft region. Given that the CMW arguments [29] 
were derived explicitly by considering the subleading behaviour of strongly ordered (coherent) parton 
showers, we do not perceive of this is as any major drawback. Instead, one should merely be aware of 
the slight shifts in the NLO corrections that result from applying it or not, recalling that a rescaling of 
A by the CMW factor ~ 1.5 is within the ordinary factor 2 variation of the renormalization scale that is 
often taken as a standard for uncertainty estimates. 

The shifts caused by CMW rescaling and/or by renormalization-scale prefactors are of course fully 
taken into account in our implementation in the VINCIA code, and are thus reabsorbed into the one-loop 
matching coefficient at the matched order, stabilizing the prediction. Differences at higher orders will 
result from the fact that the CMW rescaling, if applied, is used for all shower branchings, while the NLO 
correction derived here is only applied at the Z — )• 3 stage of the calculation. 

Because smooth p^ -ordering is the default in VINCI A we wish to understand this case as best as we 
can, and therefore we include some further comparisons with non-default settings in fig. 1 1 . 

In fig. 11a, we modify the normalization of the evolution variable from the VINCI A default = 
to the ARIADNE choice Q'^ = p\ . Though the normalization factor cancels in the Pimp factor for 
sequential gluon emissions, it is relevant for deciding the relative ordering between gluon emissions and 
gluon splittings. As this plot shows, however, the modification only produces quite small differences 
in the NLO correction factor, and with the "wrong" sign. Thus, we retain N_i_ = 4 as the default in 
VINCIA. In fig. 1 lb, we change the evolution variable for gluon splittings to be the same as that for 
gluon emissions, i.e., p^, with similar conclusions as for the previous variation. In fig. 11c, we switch 
off the Ariadne factor for gluon splittings. We notice that the NLO corrections get slightly larger. There 
is no change along the diagonal yij = yjk since the Ariadne factor is unity there, but along the edges of 
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the plots, the NLO corrections become larger, which further motivates the choice of keeping the Ariadne 
factor switched on by default in VINCIA. 

The overall result is that the infrared limits are generally well-behaved for p± evolution with ^ps oc 
p±. Remaining differences amount to small finite shifts of order 10%-20% away from unity. At that 
level, the effective finite terms of the antenna functions also play a role, hence it is too early to draw 
definite conclusions just based on the plots presented here. The impact of finite terms will be studied in 
section 5 in the context of matching to the LO matrix elements for Z — 4 partons, which effectively 
fixes the finite terms with respect to the pure-shower answers studied here. 

4.3 Tables of Infrared Limits 

The results of the preceding subsections on the infrared limits of the pole-subtracted matrix elements 
and of the Sudakov integrals generated by the various evolution-scale choices are collected here, in 
parametric form, for easy reference. The renormalization terms, Vs^, are not included. Tab. 3 expresses 
the limits of SVirtual , while tab. 4 contains the Sudakov-integral limits. 



SVirtual 


soft 




Ca + InpL 


hard collinear 


-^LCa + InpL 



Table 3: Limits of SVirtual , with L denoting ln{y), and Sqg/s = Sgq/s = y — )• 0, omitting an overall 
factor of as/27r. 





strong 


smooth 




p± 


soft 


(^L^ -iL + f^CA + lnpL 


(^L^ -iL-f^CA + lnpL 


-/3oL 


hard collinear 


-\LCa + InpL 


(^-IL-^^Ca + ^ufL 






soft 




(l' + Il--4)ca 




hard collinear 


-\LCa + IripL 


(^-lL-^)CA + lnFL 





Table 4: Limits of strong and smooth p_i_ and mr) ordering, with L denoting ln(y), with Sqg = Sgq = 
y — 0. Non divergent terms, such as vr^ have been omitted in the calculation of Vs^, and the renormal- 
ization term in V^z is set to zero. An overall factor of as/2ii is suppressed. 



5 Results including both LO and NLO corrections 

In the preceding section, we focussed on deriving the analytic forms of the shower integrals and compar- 
ing their infrared limits to the matrix-element expressions. It is now time to include also the finite terms 
arising from matching to the 4-parton tree-level matrix element, expressed by the 5 A terms in eq. (90). 
Our ultimate aim in this section is to include the full leading-colour one-loop corrections through second 
order in (i.e., up to and including Z — > 3 partons) and combine these with the full-colour tree-level 
corrections through third order in a<j (i.e., up to and including Z — )• 5 partons, the default in viNCiA). 
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However, since we shall perform the 5 A integrals numerically, adding those terms to the analytic ones 
derived in the previous section, we first wish to examine the numerical size and precision required on 
the 6 A terms themselves. 

5.1 Finite antenna terms and LO matching corrections 

Finite-term variations of the antenna functions (and in particular fixing the finite terms via unitary LO 
matching corrections, such as is done in VINCI A [18]) will affect the terms generated by the 3—7-4 
Sudakov expansions in the following way. Larger finite terms will cause an increased amount of 3 — )• 4 
branchings, which in turn will decrease the associated Sudakov factor (in the sense of driving it closer to 
zero). This will feed into the NLO correction factor, which compensates and drives the final answer back 
towards its NLO-correct value. (Note that similar variations will not occur for the 2 — 3 branching step, 
since we treat that as fixed to the LO 3-parton matrix element throughout.) This feedback mechanism is 
encoded in the 6 A terms in eq. (90). 

Following the reasoning above, we should expect larger antenna finite terms to increase the NLO 
correction factor (since, to stabilize the 3-parton exclusive rate, it must compensate for losing more 
3-parton phase-space points to 4-parton ones), and vice versa: smaller finite terms should result in a 
decrease of the NLO correction. At the pure-shower level (i.e., without LO matrix-element corrections 
to fix the finite terms), this is illustrated by fig. 12. For ease of comparison, all plots use the CMW 
rescaling of Aqcd> /Ups = P±, np = 5, and as{mz) = 0.12. The default antenna functions in VINCIA^*^ 
are shown in the middle panes, for strong (top row) and smooth (bottom row) ordering, respectively. A 
variation with smaller finite terms for the 3 — 4 antenna functions is shown to the left, and one with 
larger finite terms on the right. As expected, the NLO correction factors react by becoming lower for 
smaller finite terms and higher for larger finite terms, for both strong and smooth ordering. 

We emphasise that the plots in fig. 12 are shown purely for illustration, to give a feeling for the 
changes produced by finite-term variations. In the actual matched shower evolution, the constraint 
imposed by matching to the LO 4-parton matrix elements fixes the finite terms, via the unitary procedure 
derived in [18], which was briefly recapped in section 3.1. The effective finite terms then depend on 
the full LO 4-parton matrix elements, and have a more complicated structure than the simple antenna 
functions we have so far been playing with. We shall therefore not attempt to integrate them analytically, 
but prefer instead to let VINCIA compute a numerical MC estimate for us. 

Each point in that MC integration will involve computing at least one LO 4-parton matrix element, 
hence it is crucial to know how many points will be needed to obtain sufficient accuracy. Since ev- 
erything else is handled analytically, this will be the deciding factor in determining the speed of the 
NLO-corrected algorithm. We shall perform a speed test below in section 5.4, but first we need to 
determine the accuracy we need on the 6 A integral. 

A first analytic estimate of the size of the 6 A terms can be obtained by simply computing the ones 
produced by switching from GGG to the viNCiA default antennae (summed and averaged over helici- 

'"Note that VINCIA was recently updated with a set of hehcity-dependent antenna functions [25], so the defaults shown here 
are not identical to the GGG ones, but are instead helicity sums/averages over the functions defined in [25], 
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Figure 12: NLO correction factor for strong (top row) and smootii (bottom row) p^-ordering, for MIN 
(left), VINCIA default (middle), and MAX (right) antenna functions. We use ^ur = p± combined with 
CMW rescaling, as = 0.12, and gluon splitting in mgq. 



ties [25]), with the following 0(1) finite-term differences: 

qg ^ qgg : F^'Z'^ - F^^^, = 1.5 - (2.5 - y^, - 0.5y,-,) 
qg ^ Qq'q' : F^^-^^ - F^^,^ = 0.0 - (-0.5 + y,,) 



-l + yyij + 0.5yjk , (153) 
0.5 -y^,, (154) 



with FEmit and Fgpiit defined in eqs. (4) and (5). The 6A terms produced by these differences are 
plotted in fig. 13, for strong ordering in rriD (left) and p± (center), and for smooth ordering in p± 
(right), respectively. As expected, they do come out to be numerically subleading, roughly of order 
as/ (2vr), relative to LO (unity), yielding corrections ranging from a few permille to about a percent of 
the LO result. 

Finally, in fig. 14, we include the full LO 4-parton matrix elements and plot the distribution of 
numerically computed 6 A terms during actual VINCIA runs, for 100,000 events. The result is now 
represented by a one-dimensional histogram, with SA on the x-axis and relative rate on the y-axis. On 
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Default Settings - 
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Figure 14: Distribution of the size of the 5 A terms (normaUzed so the LO result is unity) in actual VIN- 
CIA runs. Left: linear scale, default settings. Right: logarithmic scale, with variations on the minimum 
number of MC points used for the integrations (default is 100). 



the left-hand pane, the 6A distribution with default settings is shown on a linear scale, while the right- 
hand pane shows the same result on a logarithmic scale, including variations with higher numerical 
accuracy. 

As mentioned above, the integration is done by a uniform Monte Carlo sampling of the 6 A inte- 
grands. We require a numerical precision better than 1% on the estimated size of the term (relative to 
LO) and, by default, always sample at least 100 MC points for each antenna integral. In the left-hand 
pane of fig. 14, we see that, even with the full 4-parton LO matrix-element corrections included, the size 
of the 6 A terms remains below one percent for the vast majority of 3-parton phase-space points. 

On the logarithmic scale in the right-hand pane of fig. 14, however, it is evident that there is also 
a tail of quite rare phase-space points which are associated with larger 6A corrections. Numerical 
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investigations reveal that this tail is mainly generated by the integrals over the g ^ qq terms, in particular 
in phase-space points in which the gluon is coUinear to one of the original quarks. This agrees with our 
expectation that these terms are the ones to which the pure shower gives the "worst" approximation, and 
hence they are the ones that receive the largest matrix-element corrections. As a test of the numerical 
stability of the NLO corrections for these points, we increased the minimum number of MC points 
used for the 5A integration from the default 100 (shown with "+" symbols) to 400 ("x" symbols) 
and 1600 ("*" symbols), cutting the expected statistical MC error in half with each step, at the cost of 
increased event-generation time. Though we do observe a slight broadening of the distribution between 
the default and the higher-precision settings, the shifts should be interpreted horizontally and remain 
well under the required percent-level precision with respect to LO. The default settings are therefore 
kept at a minimum of 100 MC points, though we note that future investigations, in particular of more 
complicated partonic topologies, may require developing a better understanding of, and perhaps a better 
shower approximation to, these integrals, especially the g ^ qq contribution. 

For completeness, we note that the runs used to obtain these distributions were performed using the 
new default "Nikhef" tune of viNClA's NLO-corrected shower, which will be described in more detail 
in the following subsection. Parameters for the tune are given in appendix D. 

5.2 LEP Results 

Since we have restricted our attention to massless partons in this work, we shall mainly consider the 
light-flavour-tagged event-shape and fragmentation distributions produced by the L3 experiment at LEP 
for our validations and tuning, see [55]. We consider three possible VINCIA settings: 

• New default (NLO): uses two-loop running for Ug, with CMW rescaling of Aqcd- From the 
comparisons to event-shape variables presented in this section, we settled on a value of as{Mz) = 
0.122. A few modifications to the string-fragmentation parameters were made, relative to the old 
default, to compensate for differences in the region close to the hadronization scale. The revised 
parameters are listed in appendix D, under the "Nikhef" tune. 

• New default (NLO off). Identical to the previous bullet point, but with the NLO correction factor 
switched off. 

• Old default (LO tune): uses one-loop running for Ug, without CMW rescaling of Aqcd» and 
as{Mz) = 0.139. The string-fragmentation parameters are those of the "Jeppsson 5" tune, see 
appendix D. 

The three main event-shape variables that were used to determine the value of as{Mz) are shown in 
fig. 15, with upper panes showing the distributions themselves (data and MC) and lower panes showing 
the ratios of MC/data, with one- and two-sigma uncertainties on the data shown by darker (green) and 
lighter (yellow) shaded bands, respectively. The Thrust {left) and C-parameter (middle) distributions 
both have perturbative expansions that start at 0{as) and hence they are both explicitly sensitive to the 
corrections considered in this paper. The expansion of the D parameter (right) begins at O(a^). It is 
sensitive to the NLO 3-jet corrections mainly via unitarity, since all 4-jet events begin their lives as 3-jet 
events in our framework. It also represents an important cross-check on the value extracted from the 
other two variables. 
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For a pedagogical description of the variables, see [55]. Pencil-like 2-jet configurations are to the left 
(near zero) for all three observables. This region is particularly sensitive to non-perturbative hadroniza- 
tion corrections. More spherical events, with several hard perturbative emissions, are towards the right 
(near 0.5 for Thrust and 1.0 for C and D). The maximal r = 1 — T for a 3-particle configuration 
is r = 1/3 (corresponding to the Mercedes configuration), beyond which only 4-particle (and higher) 
states can contribute. This causes a noticeable change in slope in the distribution at that point, see 
fig. 15a. The same thing happens for the C parameter at C = 3/4, in fig. 15b. The D parameter is 
sensitive to the smallest of the eigenvalues of the sphericity tensor, and is therefore zero for any purely 
planar event, causing it to be sensitive only to 4- and higher-particle configurations over its entire range. 

Both the new NLO tune (solid blue line with fiUed-dot symbols) and the old LO one (dashed magenta 
line with open-triangle symbols) reproduce all three event shapes very well. With the NLO corrections 
switched off (solid red line with open-circle symbols), the new tune produces a somewhat too soft 
spectrum, consistent with its low value of as{Mz) not being able to describe the data without the 
benefit of the NLO 3-jet corrections. 

As a further cross check, we show two further event-shape variables that were included in the L3 
study in fig. 16: the Wide and Total Jet Broadening parameters. By/ and Bt, respectively. These have 
a somewhat different and complementary sensitivity to the perturbative corrections, compared to the 
variables above, picking out mainly the transverse component of jet structure. They ai^e equal at 0{as), 
but Bt receives somewhat larger 0{a1) corrections than By/- Again, we see that both the old (LO) and 
new (NLO) defaults are able to describe the data, and that the spectrum with the new default value for 
as{Mz) is too soft if the NLO corrections are switched off. 

Finally, as an aid to constraining the Lund fragmentation-function parameters, the L3 study also 
included two infrared-sensitive observables: the charged-particle multiplicity and momentum distribu- 
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Figure 17: L3 light-flavour fragmentation observables: charged-track multiplicity and momentum dis- 
tribution. 



tions, to which we compare in fig. 17, with the momentum fraction defined as 




(155) 



There is again no noteworthy differences between the old and new default tunes. 

Having determined the value of as{Mz) and the parameters of the non-perturbative fragmentation 
function, we extended the validations to include a set of jet-rate and jet-resolution measurements by the 
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Table 5: (x^) values for: Top: L3 light-flavour event shapes and fragmentation variables [55], and LEP 
average meson and baryon fractions [61,62]. Bottom: Durham kx n-jet rates, r„j, and jet resolutions, 
y-ij, measured by the ALEPH experiment [56]. For the latter, the <^x^) calculation was restricted to the 
perturbative region, ln(y) > —8. A flat 5% theory uncertainty was included on the MC numbers. Both 
default PYTHIA and the VINCIA (LO) tune use as{mz) = 0.139 while the VINCIA (NLO) tune uses 
as{mz) = 0.122. 



ALEPH experiment [56] (now without the benefit of light-flavour tagging), using the standard Durham 
algorithm for e~^e~ collisions [57], as implemented in the FAST JET code [58]. We also compared to 
default PYTHIA 8 and, for completeness, checked that the relative production fractions of various meson 
and baryon species were indeed unchanged relative to the old VINCIA default. 

Rather than presenting all of this information in the form of many additional plots, tab. 5 instead 
provides a condensed summary of all the validations we have carried out, via (x^) values for each of 
the models with respect to each of the LEP distributions, including a flat 5% "theory uncertainty" on the 
MC numbers. Already from this simple set of values, it is clear that the LO models/tunes are already 
doing very well'^ This agreement, however, comes at the price of using a very large ("LO") value for 
Os, which is not guaranteed to be universally applicable. 

The main point of the overview in tab. 5 is that an equally good agreement can be obtained with an 
cts{mz) value that is consistent with other NLO determinations [63], specifically 

as{mz) = 0.122, (156) 

once the NLO 3-jet corrections are included. This should carry over to other NLO-corrected pro- 
cesses, and hence the fragmentation parameters we have settled on should be applicable to future NLO- 
corrected studies with VINCIA, and can also serve as a starting point for NLO-level matching studies 
with PYTHIA 8. In the latter context, the 2-loop running in particular could be retained, while the 
soft fragmentation parameters would presumably have to be somewhat readjusted to absorb differences 
between VINCIA and PYTHIA 8 near the hadronization scale'^. 

"Both VINCIA and PYTHIA are known to give quite good fits to LEP data [18,25,59,60]. For comparisons including other 
generators and tunes, see mcplots.cern.ch. 

'"The differences in soft fragmentation parameters between existing LO VINCIA and PYTHIA-8 tunes could be used as an 
initial guideline for such an effort, see, e.g., appendix D. 
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5.3 Uncertainties 



As in previous versions of VINCI A, we use the method proposed in [18] to compute a comprehensive 
set of uncertainty bands, which are provided in the form of a vector of alternative weights for each 
event. Each set is separately unitary, with average weight one'^^. The difference with respect to previous 
versions is that each variation now benefits fully from the inclusion of NLO corrections. 

When setting the parameter Vincia : uncertaintyBands = on, the uncertainty weights are 
accessible through the method 

double vincia . weight ( int i=0); 

with i = corresponding to the ordinary event sample, normally with all weights equal to unity, and the 
following variations, for i =: 

1. Default: since the user may have chosen other settings than the default, the default is included as 
the first variation. 

2. alphaS-Hi: all renormalization scales are decreased to /i = ^idei/k^, where //def = P± for gluon 
emission and ^^^i = 'niqq for gluon splitting. The default value of /c^ = 2 can be changed by 
the user. A second-order compensation for this variation is provided by the renormalization-scale 
sensitive term Vs^. 

3. alphas -Lo: /i = /ijof * k^, with similar comments as above. 

4. ant-Hi: antenna functions with large finite terms (MAX [25]). This variation is already compen- 
sated for by LO matching and is not explicitly affected by the NLO corrections. 

5. ant-Lo: antenna functions with small finite terms (MIN [25]), with similar comments as above. 

6. NLO-Hi: branching probabilities are multiplied by a factor (1 + Ug) to represent unknown (but 
finite) NLO corrections. Is canceled by NLO matching. 

7. NLO-Lo: branching probabilities are divided by a factor (1 + a^). Is canceled by NLO matching. 

8. Ord-pTl: ordering with alternative (user-definable) normalization factor in Q"^^ = N_i_p\. 
Compensated at first order by LO matching, and at second order (Sudakov corrections) by NLO 
matching via ordering-sensitive terms in Vjj^. 

9. Ord-mD: smooth mD ordering, with similar comments as above. 

10. NLC-Hi: qg emission antennae use Ca as color factor. Compensated at first order by LO match- 
ing. Not affected by NLO matching since we only operate at leading color. 

11. NLC-Lo: qg emission antennae use 2Cf as color factor, with similar comments as above. 

"vincia currently does not attempt to give a separate estimate of the uncertainty on the total inclusive cross section. The 
uncertainties it computes only pertain to shapes of distributions and the effects of cuts on the total inclusive rate. 
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We emphasize that these variations are not all independent (for instance the and NLO variations are 
highly correlated) and hence the corresponding uncertainties should not be summed in quadrature. In the 
VINCIAROOT plotting tool included with VINCIA, the uncertainty band is constructed by taking the max 
and min of the variations. See the viNCiA HTML manual for more information about the uncertainty 
bands and [18] for details on their algorithmic construction. 

5.4 Speed 

Although the CPU time required by matrix-element and shower/hadronization generators is still typi- 
cally small in comparison to that of, say, full detector simulations, their speed and efficiency are still 
decisive for all generator-level studies, including tuning and validation, parameter scans, development 
work, phenomenology studies, comparisons to measurements corrected to the hadron level, and even 
studies interfaced to fast detector simulations. For this wide range of applications, the high-energy sim- 
ulation itself constitutes the main part of the calculation. An important benchmark relevant to practical 
work is for instance whether the calculation can be performed easily on a single machine or not. 

Higher matched orders are characterized by increasing complexity and decreasing unweighting ef- 
ficiencies, resulting in an extremely rapid growth in CPU requirements (see e.g. [25]). At NLO, the 
additional issues of negative weights and/or so-called counter-events can contribute further to the de- 
mands on computing power. With this in mind, high efficiencies and fast algorithmic structures were a 
primary concern in the development of the formalism for leading-order matrix-element corrections in 
VINCIA [18,24,25], and this emphasis carries through to the present work. We can make the following 
remarks. 

• The only fixed-order phase-space generator is the Bom-level one. All higher-multiplicity phase- 
space points are generated by (trial) showers off the lower-multiplicity ones. This essentially 
produces a very fast importance-sampling of phase-space that automatically reproduces the dom- 
inant QCD structures. 

• Likewise, the only cross-section estimate that needs to be precomputed at initialization is the total 
inclusive one. Thus, initialization times remain at fractions of a second regardless of the matching 
order. 

• The matrix-element corrected algorithm works just like an ordinary parton shower, with modified 
(corrected) splitting kernels. In particular, all produced events have the same weights, and no 
additional unweighting step is required. 

• Since the corrections are performed multiplicatively, in the form of (1 + correction), with 1 
being the LO answer, there are no negative-weight events and no counter-events. The only ex- 
ception would be if the correction becomes larger than the LO answer, and negative. This would 
correspond to a point with a divergent fixed-order expansion, in which case the use of NLO cor- 
rections would be pointless anyway. Moreover, as demonstrated by the plots in the previous sec- 
tions, our definitions of the corrections are analytically stable (and numerically subleading with 
respect to LO) over all of phase space, including the soft and coUinear regions, for reasonable 
renormalization- and evolution-variable choices. 
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LO level NLO level Time / Event Speed relative to PYTHIA 

Z — )■ Z — )• [milliseconds] / PYTHIA 8 

PYTHIA 8 2,3 2 0.4 1 

VINCIA (NLO off) 2,3,4,5 2 2.2 ~ 1/5 

VINCIA (NLO on) 2,3,4,5 2,3 3.0 ~ 1/7 



Table 6: Event-generation time in VINCIA 1.0.30 compared to PYTHIA 8.176. 



• The parameter variations described in section 5.3 can be performed together with the matching 
corrections to provide a set of uncertainty bands in which each variation benefits from the full 
corrections up to the matched orders. These are provided in the form of a vector of alternative 
weights for each event [18], at a cost in CPU time which is only a fraction of that of a comparable 
number of independent runs. 

These attributes, in combination with helicity dependence in the case of the leading-order formal- 
ism [25], allow VINCIA to run comfortably on a single machine even with full-fledged matching and 
uncertainty variations switched on. 

The inclusion of NLO corrections will necessarily slow down the calculation. The relative increase 
in running time relative to PYTHIA 8, is given in tab. 6, including the default level of tree-level match- 
ing, with and without the NLO 3-jet correction Without it (but still including the default tree-level 
corrections which go up to Z — )■ 5 partons), VINCIA is 5 times slower than PYTHIA. With the NLO 3-jet 
correction switched on, this increases only slightly, to a factor 7. For a fully showered and hadronized 
calculation which includes second-order virtual and third-order tree-level corrections, we consider that 
to still be acceptably fast. Importantly, an event-generation time of a few milliseconds per event implies 
that serious studies can still be performed on an ordinary laptop computer. 

6 Outlook and Conclusions 

In this work, we have investigated the expansion of a Markov-chain QCD shower algorithm to second or- 
der in the strong coupling, for e+e~ — 3 partons, and made systematic comparisons to matrix-element 
results obtained at the same order. Using these results, we have subjected the subleading properties 
of shower algorithms with different evolution/ordering variables and different renormalization-scale 
choices to a rigorous examination. At the analytical level, we have compared the logarithmic struc- 
tures at the edge of phase space, and at the numerical level we have illustrated the difference between 
the expanded shower algorithm and the one-loop matrix element. 

We find that the choice of -ordering, with a renormalization scale proportional to p± yields the 
best agreement with the one-loop matrix element, over all of phase space. This elaborates on, and is con- 
sistent with, earlier findings [28,29]. Using the antenna invariant mass, rriD, for the evolution variable 
still gives reasonable results in the hard regions of phase space, but leads to logarithmically divergent 
corrections for soft emissions, the exact form of which depends on the choice of renormalization vari- 

'''The numbers include both showering and hadronziation and were obtained on a single 3.05 GHz CPU, with gcc 4.7 -02, 
using default settings for PYTHIA 8 and the "Nikhef" NLO tune for VINCIA. 
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able. In the VINCIA code, we retain the option of using m£) mainly as a way of providing a conservative 
uncertainty estimate. 

With the NLO 3-jet corrections included as multiplicative corrections to the shower branching prob- 
abilities, we find that we can obtain good agreement with a large set of LEP event-shape, fragmentation, 
and jet-rate observables with a value of the strong coupling constant of as{Mz) = 0.122. This is in 
strong contrast with earlier (LO) tunes of both PYTHIA and VINCI A which employed much larger values 
~ 0.14 to obtain agreement with the LEP measurements. 

This paper is intended as a first step towards a systematic embedding of one-loop amplitudes within 
the VINCIA shower and matching formalism. To arrive at a full-fledged prescription, this will need to 
be extended to hadron collisions, ideally in a way that allows for convenient automation. A first step 
towards developing the underlying shower formalism for pp collisions was recently taken [64], and more 
work is in progress. 

In addition, further studies should be undertaken of the impact of unordered sequences of radia- 
tion that can occur for the smooth-ordering case (it may be necessary to adopt a strategy similar to 
the truncated showers of the MC@NLO approach), and the mutually related issues of total normaliza- 
tion and how much of the (hard) corrections are exponentiated (similar to the differences between the 
POWHEG and MC @ NLO formalisms, but here occurring at one additional order, where the relevant total 
normalization is the NNLO one). Finally, it would be interesting to develop an extension of this formal- 
ism that would allow second-order-corrected antenna functions to be used at every stage in the shower, 
thereby upgrading the precision of the all-orders resummation, a project that would involve examining 
the second-order corrections to branchings of qg and gg mother antennae as well. We look forward to 
following up on these issues in the near future. 
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A Infrared singular operators 

Here we list the IR singularity operators from [16,27,51] as they are used in section 3. 



-^qq ) A* / ^qq) 
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e^T r 1 5 
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B One-Loop Amplitudes 



B.l Renormalization 

Because a detailed derivation of the calculation of Z — )• 3 jets can be found in [26] we restrict ourselves 
to listing the result in form that is convenient for our purpose. Divergences are regulated using dimen- 
sional regularization with d = 4 — 2e. Our results, before ultraviolet renormalization, are cross-checked 
with [26] where one must undo the renormalization in their case. In order to cancel the ultraviolet poles 
we need to renormalize the coupling according to (see also section 3.4) 
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where 



/3o 



lliVc - 2nF 
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and Se = (47r)'^ exp(— e7£;) contains the factors characterizing the MS scheme. Due to the renormal- 
ization, the leading order calculation will generate a term quadratic in as(/ufj). 
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which directly cancels the ultraviolet poles of the next-to-leading order calculation. 



B.2 One-loop Matrix Element 

The fixed-order expression relevant to matching in the VINCIA context is the one-loop matrix element 
normalized by the tree-level one. We decompose this into leading- and subleading-colour pieces, as 
follows: 

2 Re 



' {LC + QL + SL) , (163) 



|M0|2 2tt 

with the LC piece containing the Ca part of the gluon loops, the QL one containing the quark loops, 
proportional to npTfj, and the SL piece containing the subleading gluon-loop corrections, proportional 
to —1/Nc- As usual in MC applications, we usually refer to "Leading Colour" as including both the 
Nc and Tr pieces. These are both associated with so-called planar colour flows that are simple to 
relate to the colour-flow representations used in Monte Carlo event generators, see e.g. [1,65]. The 
subleading-colour piece is included below for completeness, but has so far been left out of the NLO 
matching corrections implemented in the viNCiA code. 

The notation of the infrared pole structure of these terms has been written similar to the integrated 
antenna in [16], with the difference that we have chosen to expand the scale factor fi in the integrated 
antenna terms in order to obtain explicitly dimensionless logarithms. 
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The quark has been labelled 1, the anti-quark 2 and the gluon 3. 
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and the infrared singular operators, given in appendix A. 

With the one-loop matrix element expressed in this form, cancellation of the infrared poles against 
the integrated antennae (see below) coming from the shower will be particularly simple and will yield 
an expression purely dependent on the renormalization scale, fiR, and on the kinematic invariants S12, 
S23, and S12, but not on the scale factor /i. 
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C Antenna integrals 



In this appendix we list the results of antenna integrals over phase space corresponding to the various 
evolution variables. 

C.l Strong Ordering Gluon Emission 

The expressions for a gluon emitting antenna is given in eq. (4). With a redefinition the same antenna 
function reads 
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where yi, 7/2 correspond to yij, of eq. (6), respectively. Recall that the last three terms serve to give 
a flexible and explicit way of tracking extra non-singular terms in antennae. The phase space integral 
over these antenna, as determined by the evolution variable, can be written as 
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All antenna integrals in eq. (68) have been written in such a way that they are integrated over their whole 
invariant mass plus a correction term running from the evolution variable to the total invariant mass. The 
integrals running over the whole invariant mass contain singular regions and therefore poles while the 
correction terms yield finite corrections. These finite corrections are discussed per evolution variable 
below. We define the integrals 
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h = 2?Qc [Coo + CoiyiiQ^ C) + Cioy2(Q', 0] • 
So that, in these terms, the results read 



1 /■'"/if 1 / ^ \ 



where 



Ki = l, K2 = -2 Ks = 2, Ki = -6jg-5Kg, = I. 
We now turn to specific cases. 
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C.1.1 Dipole Virtuality 

The results for the individual contributing parts of the antenna function as defined in eq. (172) - eq. (176) 
withe = ^.'"'"^'.^"/^^'^ . are 
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In the case of integration over the 3 — )• 4 splittings, the definition of the integrals remains the same. 
Only the definition of g changes with 
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C.1.2 Transverse Momentum 

The results for the individual contributing parts of the antenna function as defined in eq. (172) - eq. (176) 
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With ?/i = and = s. In the case of the 3 — )• 4 splittings the only adaptation takes place in the 

former definition where is set equal to Sqg or Sgq dependent on which dipole is being integrated 
over. 



C.1.3 Energy ordering 

The results for this evolution parameter are 
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with A used as a cut-off on 4p]_ and y3 = i^aai^Ml!. 
C.2 Strong Ordering Gluon Splitting 

The branching of a gluon splitting into a quark antiquark pair can only take place at the 3 — 4 level 
splitting. The generation of a gluon splitting takes place through an alternative form of phase space 
generation than the discussed mo, p± and En variables. Instead phase space is sampled in a triangular 
surface comparable to ordering, yet in this case using only one cutoff, the generated at the 2 — > 3 
level, to avoid the singular region of the gluon sphtting antenna. The gluon spUtting antenna is given by 

aq/ggiyi,y2) = 4- f ^^^^ + - + Coo + Com + Cioy2) • (180) 

Because the integration surface is similar for all evolution types only depending on the cutoff the 
integration is demonstrated for all types 
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VE (-6Cio(3 + (-3 + yE)yE + yE{-9 + ^yE + 6Coi(-3 + 2yE))) - 12 In(yE)] 

Where the factor a half has been added for the sake of consistency with respect to the treatment of gluon 
emission. The factor needs to be replaced by either Sqg or Sgq dependent on which dipole is being 
integrated. 
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C.3 Smooth Ordering Gluon Emission 



The phase space integral in the case of smooth ordering differs from strong ordering by allowing inte- 
gration over the whole phase space region. The inclusion of a damping factor regulates the accessible 
region of phase space which generates a different phase space occupancy than in the case of strong 
ordering. A general form for smooth ordering integration of a gluon emission antenna is 
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Where we use the definition of eq. (169) with Si = yirnjj^, denotes the branching scale and Qe 
indicates the evolution variable used for gluon emission. We define the following integrals 
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So that, in these terms, the results read 
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where 

Ki = l, K2 = -2 K3 = 2, K4 = -6ig-6Kg, = I. 
We now turn to specific cases. 
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C.3.1 Smooth mass ordering 



The only term from eq. (182) that requires specification is the damping factor 
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The computation of the individual antenna parts will require separating the phase space triangle in two 
regions (si > S2 and vice versa) in order to make the damping factor definite. After summing over these 
two regions we obtain the following values for gluon emission contributions 
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C.3.2 Smooth transverse momentum ordering 



In the case of smooth ordering for transverse momentum we find the following result for the ordering 
requirement 



(192) 



Where rnjj^ should be replaced by Sqg or Sgq dependent on the dipole of integration. In combination 
with eq. (182) we find the following results for the partial gluon emission antenna parts 
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C.4 Smooth Ordering Gluon Splitting 

Additionally we also need to consider the gluon splitting antenna function for smooth ordering. Similar 
to the strong ordering case, the separate generation of gluon splitting variables allows for a new choice 
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for evolution variable and thereby a different phase space surface. As in the case of gluon emission we 
allow for integration over the whole phase space, using the damping factor to limit the accessible area. 
A general notation is the following 

"^""W ds2^^^^f^a5/gg(si,S2). (193) 

^IK Jq% Jo + ^3 

With the definition for the gluon splitting antenna as in eq. (180). 
C.4.1 Dipole mass 

With the gluon splitting antenna as defined in eq. (180) and the phase space integral eq. (193) we find 
the following result for the dipole mass evolution type 

G =^ (-12(4 + 3(-l + 2C7oo + Cio)yi4 + 3(-l + 2Coo - 2Coi + 2Cio + 2Cio)y|sp 



72sj, 



+ (-2 - 6Coi + 3Cio)yf ) arctanh 



Sp 



2yl + sp 



+ sp ((-13 + 18Coo + 6C7oi + 6Cio)4 + 3(-4 + 12Coo - 6Coi + 9Cio)2/|sp 
-6(2 + 6Coi - 3Cio)y| + 184 In \ \ (194) 



^3 



With y3 = ^ '!am{sqg,sgq) Sp = max(sgg, Sgq)/ s. Note that the gluon splitting antenna has been 
defined with the singularity in 2/2 which determines the form of the damping factor. 

C.4.2 Transverse momentum 

With the gluon splitting antenna defined in eq. (180) and the phase space integral eq. (193) we find the 
following result for the transverse momentum evolution type 

G =^ (4 (8 (-13 + 6C01 - 32/i) + 72Coo(2 + yf) 



3yi(2yi + 6C7oi(2 + yl) - 3Cio)(6 + yf)) + 72y| arctanh 



2 



.2 + yi, 

+ 3 (128 + yl (-24Coo(4 + yf) - 3C7io(4 + ylf + 2yi(6 + yf + 3Coi(4 + y|)))) 
xlogfi±?^)) (195) 



Withy? = 4^^229£££. 
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D NLO l\ine Parameters 



In tab. 7 below, we list the perturbative and non-perturbative fragmentation parameters for the Nikhef 
NLO tune of VINCIA. For reference, we compare them to the current (LO) default Jeppsson 5 tune, 
which was used for comparisons to LO viNCiA in this paper. 





NLO Tune 


LO Tune 




Parameter 


(Nikhef) 


(Jeppsson 5) 


Comment 


! * alphas 








Vincia : alphaS value 


= 0.122 


= 0.139 


! alphas (mZ) value 


Vincia:alphaSkMu 


= 1.0 


= 1.0 


! Renormalization-scale prefactor 


Vincia:alphaSorder 


= 2 


= 1 


! Running order 


Vincia:alphaSmode 


= 1 


= 1 


! muR = pT:emit and Q: split 


Vincia:alphaScmw 


= on 


= off 


! CMW rescaling of Lambda on/off 



! * Shower evolution and IR cutoff 






Vincia : evolutionType 


= 1 


= 1 


! pT-evolution 


VinciaiorderingMode 


= 2 


= 2 


! Smooth ordering 


Vincia:pTnormalization 


= 4. 


= 4. 


! QT = 2pT 


Vincia : cutoffType 


= 1 


= 1 


! Cutoff taken in pT 


VinciaxutoffScale 


= 0.8 


= 0.6 


! Cutoff value (inGeV) 



! * Longitudinal string fragmentation parameters 



StringZ:aLund 


= 0.40 


= 0.38 


! Lund FF a (hard fragmentation supp) 


StringZ:bLund 


= 0.85 


= 0.90 


! Lund FF b (soft fragmentation supp) 


StringZ: aExtraDiquark 


= 1.0 


= 1.0 


! Extra a to suppress hard baryons 



! * pT in string breakups 



String 


'PTisigma 


= 0.29 


= 0.275 


String 


;PT:enhancedFraction 


= 0.01 


= 0.01 


String 


;PT: enhancedWidth 


= 2.0 


= 2.0 


! * String breakup flavour parameters 




String 


;Flav:probStoUD 


= 0.215 


= 0.215 


String 


iFlav:mesonUDvector 


= 0.45 


= 0.45 


String 


;Flav : mesonS vector 


= 0.65 


= 0.65 


String 


;Flav:mesonC vector 


= 0.80 


= 0.80 


String 


;Flav:probQQtoQ 


= 0.083 


= 0.083 


String 


;Flav:probSQtoQQ 


= 1.00 


= 1.00 


String 


;Flav:probQQltoQQO 


= 0.031 


= 0.031 


String 


;Flav:etaSup 


= 0.68 


= 0.68 


String 


;Flav : etaPrimeS up 


= 0.11 


= 0.11 


String 


'Flav : decupletS up 


= 1.0 


= 1.0 



Table 7: Parameters of the "Nikhef" NLO tune. 



Soft pT in string breaks (in GeV) 
Fraction of breakups with enhanced pT 
Enhancement factor 



! Strangeness-to-UD ratio 

! Light-flavour vector suppression 

! Strange vector suppression 

! Charm vector suppression 

! Diquark rate (for baryon production) 

! Optional Strange diquark suppression 

! Vector diquark suppression 

! Eta suppression 

! Eta' suppression 

! Optional Spin-3/2 Baryon Suppression 
compared to those of the "Jeppsson 5" LO tune. 
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